fcT'L  INST.  OF   STAND  &  TEf"  f.:]f. 


IIIDS    blfiT17 


NIST 
PUBLICATIONS 


Nisr 


United  States  Department  of  Commerce 

Technology  Administration 

National  Institute  of  Standards  and  Technology 


N/ST  Technical  Note  1385 


NIST  Technical  Note  1385 

Global  Position  System  Receivers 
and  Relativity 


Neil  Ashby 
Marc  Weiss 


Time  and  Frequency  Division 

Physics  Laboratory 

National  Institute  of  Standards  and  Technology 

325  Broadway 

Boulder,  Colorado  80303-3328 


March  1999 


'4rES 


U.S.  DEPARTMENT  OF  COMMERCE,  William  M.  Daley,  Secretary 
TECHNOLOGY  ADMINISTRATION,  Gary  R.  Bachula,  Acting  Under  Secretary  for  Technology 
NATIONAL  INSTITUTE  OF  STANDARDS  AND  TECHNOLOGY,  Raymond  G.  Kamnier,  Director 


National  Institute  of  Standards  and  Technology  Technical  Note 
Natl.  Inst.  Stand.  Technol.,  Tech.  Note  1385,  52  pages  (March  1999) 

CODEN:NTNOEF 


U.S.  GOVERNMENT  PRINTING  OFFICE 
WASHINGTON:  1999 


For  sale  by  the  Superintendent  of  Documents,  U.S.  Government  Printing  Office,  Washington,  DC  20402-9325 


Contents 

1 .  Introduction 1 

2.  Theory 2 

2.1  Pseudorange   2 

2.2  Phase  Time  2 

2.3  Relevant  Relativity 3 

2.4  Common  Misunderstandings 4 

2.5  User  Corrections 6 

2.6  Relativistic  Doppler  Effect  and  Use  of  the  GPS  Carrier 7 

2.7  Rotation  Matrix 8 

3.  Time  Tagging  at  Receiver 9 

3.1  Time  Tagging  at  Receiver — General  Prescription 9 

3.2  Examples  11 

Example  3.21  Multi- Channel  Receiver,  Time  Tagging  at  Receiver   12 

Example  3.211  Erroneous  Use  of  ECEF  Coordinates   16 

Example  3.212  Iteration  With  a  Succession  of  Inertial  Frames 16 

Example  3.213  The  Benefit  of  Better  Initialization 17 

4.  Time  Tagging  at  Transmitter 18 

4.1  Code  Boundary  Measurements 18 

4.2  Phase  Time 19 

4.3  Accounting  for  Receiver  Motion 20 

4.4  Time  Tagging  at  Transmitter — General  Prescription 22 

4.5  Examples  24 

Example  4.51  Time  Tagging  at  Transmitter 24 

Example  4.52  Time  Tagging  at  Transmitter,  Neglecting  Doppler  Corrections    29 

5.  Time  Transfer  Receiver 29 

Example  5.1  Earth-Fixed  Single  Channel  Timing  Receiver 32 

6.  References 32 

Appendix  A.  The  Eccentricity  Correction 35 

Appendix  B.  Linearization  Algorithm  for  Time  Tagging  at  Multichannel  Receiver 37 

Appendix  C.  Linearization  Algorithm  for  Time  Tagging  at  Transmitters 39 

Appendix  D.  Equivalence  of  Integrated  Doppler  and  One- Way  Light  Time 42 

Appendix  E.  Numerical  Check  on  First-Order  Doppler  Correction 45 


111 


GLOBAL  POSITIONING  SYSTEM  RECEIVERS  AND  RELATIVITY 


Marc  Weiss 

National  Institute  of  Standards  and  Technology 

Time  and  Frequency  Division 

325  Broadway,  Boulder,  CO  80303-3324 


Neil  Ashby 

Department  of  Physics,  Campus  Box  390 

University  of  Colorado 

Boulder,  CO  80309-0390 


ABSTRACT 

We  illustrate  the  general  methods  for  applying  relativistic  corrections  needed  by  a  GPS 
receiver  in  providing  time  or  position  to  a  user.  We  focus  on  estimating  the  geometric  range 
delay,  the  time  GPS  signals  require  to  propagate  from  the  transmitter  to  the  receiver  in 
vacuum.  We  discuss  several  common  uses  of  GPS  including  time  tagging  at  the  receiver, 
time  tagging  at  the  transmitter,  and  earth-fixed  timing  receivers.  We  include  examples 
which  illustrate  these  cases  numerically. 
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1.  INTRODUCTION 

A  GPS  receiver  must  make  two  corrections  that  are  related  to  relativity  in  order  to 
provide  time  or  position  to  a  user.  We  discuss  these  corrections  and  focus  mostly  on 
estimating  the  geometric  range  delay  to,  the  time  for  GPS  signals  to  propagate  from  the 
transmitter  to  the  receiver  in  vacuum.  Proper  estimation  of  t^  is  essential  for  solving  for 
position  or  time.  This  is  an  application  of  the  relativistic  principle  of  the  constancy  of  the 
speed  of  light,  which  states  that  electromagnetic  signals  travel  in  Euclidean  straight  lines 
with  speed  c  relative  to  an  inertial  reference  frame.  We  present  a  few  cases  which  apply 
to  many  common  uses  of  GPS  and  illustrate  them  numerically. 

We  present  the  theory  behind  corrections  with  references  given  for  any  derivations  not 
done  here.  Through  the  derivations,  we  show  that  tlu>  Interfac-e  Control  Document  (ICD- 
GPS-200)  specifications,  as  issued  by  the  Joint  Program  Office  of  the  Global  Positioning 


System  [1],  consistently  cover  the  requirements  of  relativity  down  to  within  1  ns  or  less. 
We  show  that  the  ICD  specifications  include  relativity  corrections  with  enough  accuracy 
for  most  applications  [2,3].  In  particular,  we  discuss  the  relativistic  Doppler  effect,  the 
formula  for  its  instantaneous  magnitude,  and  its  relationship  with  typical  GPS  receiver 
operation.  We  also  address  the  use  of  carrier-phase  measurements,  which  is  not  discussed 
in  the  ICD. 

Relativistic  effects  which  are  not  usually  modelled  in  GPS  are:  (a)  a  small  effect  on 
satellite  clock  rates  due  to  the  earth's  oblateness  and  (b)  a  small  time  delay  due  to  the 
slowing  of  electromagnetic  signal  propagation  in  the  earth's  gravitational  field.  The  earth's 
oblateness  contributes  a  very  small  constant  rate  correction  for  satellite  clocks,  which  is  so 
small  (less  than  50  ps  per  day)  that  it  can  be  neglected.  It  also  causes  a  periodic  variation  in 
the  GPS  satellite  clock  time  having  twice  the  orbital  period  and  a  peak-to-peak  amplitude 
less  than  200  ps.  For  a  GPS  satellite-to-user  Hnk  the  gravitational  time  delay  is  less  than 
200  ps  [2]. 

Our  goal  is  to  give  simple  recipes,  so  that  GPS  users  may  understand  how  to  im- 
plement receiver  designs  consistently  with  the  requirements  of  relativity.  We  then  back 
up  these  prescriptions  with  full  explanations  and  derivations,  so  that  those  interested  can 
understand  where  the  prescriptions  come  from  and  what  approximations  are  involved.  The 
result  is  that  this  document  provides  a  brief  tutorial  on  how  relativity  is  apphed  for  users 
of  the  GPS  broadcast  ephemerides. 

2.  THEORY 

2.1  Pseudorange 

Generally,  a  GPS  navigation  receiver  measures  the  arrival  times,  on  a  local  clock,  of 
signals  from  four  or  more  different  satellites.  The  receiver  then  uses  these  measurements 
along  with  information  in  the  navigation  messages  to  solve  for  four  unknowns:  user  position 
{x,y,z)  and  the  receiver  clock  offset  from  GPS  time.  The  signals  from  the  sateUites  can 
be  thought  of  effectively  as  continuous  timing  signals  arriving  at  the  receiver  from  the 
satellite  clocks.  The  receiver  has  its  own  local  clock  for  comparison.  The  user  measures 
either  (1)  the  times  in  the  received  timing  signals  at  a  specific  local  clock  time,  or  (2)  the 
arrival  times  on  the  local  clock  of  a  specific  time  tag  in  the  received  timing  signals.  The 
difference  between  the  reception  time  (according  to  the  local  clock)  and  the  transmission 
time  at  one  satellite  (according  to  GPS  time),  multiplied  by  the  speed  of  fight,  is  called 
the  pseudoTU'fige.  These  pseudoranges  are  used  to  solve  for  user  position  and  time. 

2.2  Phase  Time 

In  addition  to  times  measured  on  real  clocks,  we  can  think  of  a  particular  phase  point 
in  the  wave  train  between  transmitter  and  receiver  as  nvnnerically  equal  to  some  time  r, 


which  follows  the  wave  as  it  propagates.  This  "phase  time"  r  is  transmitted  from  the 
satellite  at  SV  clock  time  tsv  and  propagates  to  the  receiver  where  it  arrives  at  receiver 
clock  time  t^.  Thus  the  pseudorange  is  c{tpi  —  tsv)- 

2.3  Relevant  Relativity 

Several  relativistic  effects  have  already  been  incorporated  into  the  GPS  system  so, 
for  the  ordinary  user  of  broadcast  ephemerides,  only  two  relativistic  corrections  must  be 
considered.  First,  the  receiver  must  apply  a  correction  to  the  transmitted  time  to  account 
for  relativistic  effects  arising  from  orbit  eccentricity  of  the  transmitting  satellite.  This  is 
the  Atr  term  defined  in  the  ICD.  Second,  the  finite  and  universally  constant  speed  c  of 
signals  propagating  in  a  vacuum  from  a  transmitter  to  a  receiver,  relative  to  an  inertial 
frame  (the  geometric  path  delay),  must  be  accounted  for.  The  receiver  must  also  account 
for  ionospheric  and  tropospheric  delay  corrections,  which  we  do  not  consider  here. 

Three  relativistic  effects  are  germane  to  GPS.  Rates  of  clocks  in  GPS  are  adjusted 
(as  they  are  for  International  Atomic  Time)  to  match  the  rate  that  clocks  would  run  on 
the  geoid  of  the  earth.  The  geoid  is  a  surface  of  constant  gravitational  potential  in  the 
rotating  frame  in  which  the  effects  (2)  and  (3)  described  below  add  to  a  constant  value. 
The  three  relativity  effects  are  as  follows. 

(1)  Constancy  of  the  speed  of  light  and  relativity  of  synchronization.  GPS  time  is 
defined  using  the  principle  of  the  constancy  of  c  to  synchronize  an  imagined  system  of 
clocks  everywhere  in  space  in  the  neighborhood  of  the  earth  (this  is  called  Einstein  syn- 
chronization). GPS  satellite  clocks  are  in  principle  adjusted  to  agree  with  this  imagined 
system  of  clocks.  This  network  of  synchronized  GPS  clocks  realizes  a  coordinate  time,  a 
system  of  self-consistent  time  markers  with  which  to  label  events.  This  definition  of  GPS 
time  requires  a  locally  inertial  coordinate  system.  GPS  time  is  thus  defined  relative  to  an 
earth-centered  inertial  coordinate  system  (an  ECIF),  but  the  rate  is  set  to  match  the  rate 
at  which  clocks  would  run  on  the  geoid.  An  ECIF  is  also  used  to  simplify  the  paths  of 
signals  propagating  from  satellites,  since,  with  sufficient  accuracy  for  GPS,  light  travels  in 
Euclidean  straight  lines  at  the  speed  c  in  vacuum  relative  to  such  inertial  frames  [2]. 

(2)  Second-Older  Doppler  shift.  A  clock  moving  with  respect  to  an  ECIF  runs  slower 
relative  to  coordinate  time  in  that  ECIF  than  if  it  were  at  rest  in  the  ECIF.  This  is  the 
time  dilation  effect  due  to  the  magnitude  of  the  relative  velocity,  sometimes  called  the 
second-order  Doppler  effect.  For  satellites  in  GPS  orbits,  the  fractional  frequency  offset 
needed  to  compensate  for  this  is  approximately  -1-8.3  x  10~^^  relative  to  the  rate  of  clocks 
on  the  earth's  geoid. 

(3)  Gravitational  frequency  shift.  A  clock  at  rest  in  a  lower  gravitational  potcMitial 
runs  slower  relative  to  coordinate  time  than  if  it  were  at  rest  in  a  higher  potential.  This 
is  called  the  gravitational  red  shift.  Thus,  standard  clocks  closer  to  the  earth  nin  slowt^- 
than  standard  clocks  farther  away,  since  the  gravitational  potential  heconu^s  nu^re  negative 


closer  to  the  earth.  Clocks  on  GPS  satellites  run  faster  than  clocks  at  rest  on  the  earth's 
surface.  Thus  GPS  satellite  clock  frequencies  need  to  be  adjusted  by  a  fraction  of  about 
—5.3  X  10~^°  relative  to  the  earth's  geoid,  to  compensate  for  this  effect. 

To  compensate  for  the  relativistic  effects  described  in  paragraphs  (2)  and  (3)  above 
for  circular  orbits,  and  as  a  consequence  of  the  requirement  that  GPS  satellite  clocks  run 
at  the  rate  that  a  standard  clock  on  the  geoid  would  run,  atomic  clocks  in  GPS  satellites 
are  given  a  fixed  fractional  frequency  offset  of  —4.4645  x  10~^^.  The  user  does  not  have 
to  be  concerned  about  these  rate  corrections. 

These  three  relativistic  effects,  however,  explain  why  the  user  must  still  apply  two 
corrections.  The  first  relativistic  correction  is  the  Atr  term  defined  in  the  ICD.  This  term 
corrects  the  satellite  vehicle  (SV)  clock  offset  due  to  any  eccentricity  in  the  vehicle  orbit. 
Eccentricity  produces  periodic  time  excursions  from  the  nominal  fixed  fractional  rate  offset 
of  —4.4645  X  10~^°.  The  second  correction  applies  to  users'  estimates  of  the  geometric  range 
delay,  the  delay  from  the  transmitter  to  the  receiver  if  the  signal  traveled  in  vacuum.  Since 
an  ECEF  frame  rotates  relative  to  any  inertial  frame,  light  does  not  travel  in  a  Euclidean 
straight  line  in  an  ECEF  frame.  It  is  incorrect  to  calculate  the  geometric  range  in  an 
ECEF  frame  using  the  magnitude  of  the  vector  difference  between  the  receiver's  position 
at  reception  and  the  transmitter  position  at  transmission.  The  distance  is  most  easily 
calculated  in  an  ECIF  where  signals  do  travel  in  Euclidean  straight  lines  at  the  speed  of 
light,  the  constant  c. 

Many  ECI  coordinate  systems,  differing  by  constant  spatial  rotations  from  each  other, 
may  serve  these  purposes.  All  frames  that  have  the  same  origin  at  the  earth's  center  and 
do  not  rotate  with  respect  to  the  "fixed"  stars  will  define  simultaneity  in  the  same  way. 
All  such  frames  are  equivalent  for  determining  the  propagation  delay.  Yet  users  usually 
need  to  refer  their  positions  to  the  earth.  For  example  they  may  need  their  coordinates  in 
the  WGS-84  coordinate  system,  the  ECEF  frame  upon  which  the  broadcast  ephemerides 
are  based.  Users  fixed  on  the  earth  who  want  time  from  GPS  often  know  their  coordinates 
as  constants  in  the  ECEF  frame.  In  fact,  most  users  would  hke  to  know  their  coordinates 
in  the  ECEF  frame. 

2.4  Common  Misunderstandings 

There  has  been  confusion  in  the  GPS  community  over  a  number  of  issues  related  to 
relativity  [5,6].  For  example,  arguments  have  been  presented  that  we  need  to  define  a 
coordinate  system  in  the  rest  frame  of  the  receiver  or  perhaps  in  that  of  the  transmitter. 
Since  measurements  are  made  in  the  receiver,  it  might  at  first  appear  necessary  to  use  a 
coordinate  system  comoving  with  the  receiver.  However,  such  arguments  overlook  the  fact 
that  relativistic  effects — indeed,  any  physical  effects — can  be  described  in  any  coordinate 
system. 

Another  argument  that  has  been  advanced  is  that  we  need  to  compute  corrections  in 


the  rotating  coordinate  system  of  the  ECEF  frame,  since  that  is  the  coordinate  system  in 
which  the  user  ultimately  wishes  to  find  his/her  position.  The  relativity  of  simultaneity 
has  sometimes  been  presented  as  a  principle  necessitating  such  constructions. 

Although  measurements  are  made  with  a  particular  instrument  and  an  instantaneous 
rest  frame  is  associated  with  that  instrument,  it  is  not  necessary  to  use  that  rest  frame 
to  analyze  the  data.  This  is  a  common  misunderstanding  underlying  many  of  the  errors 
in  the  literature:  the  belief  that  we  must  always  use  reference  frames  which  are  comoving 
with  the  transmitter  or  receiver.  The  user  does  not  need  a  coordinate  system  whose  origin 
coincides  with  the  user's  trajectory  in  spacetime  to  analyze  effects  local  to  the  user. 

It  is  a  common  practice  in  physics  to  choose  a  convenient  coordinate  system  to  solve  a 
problem.  This  is  particularly  germane  to  relativity,  where  the  principle  of  covariance  states 
that,  if  a  problem  is  understood  in  one  coordinate  system,  it  can  be  understood  in  any 
other  by  making  appropriate  transformations.  Thus,  even  though  neither  the  user  nor  the 
transmitter  are  at  rest  in  an  ECIF,  the  user  may  choose  an  ECIF  to  simplify  computations 
of  position,  then  transform  those  coordinates  to  the  ECEF  frame,  the  WGS-84  system. 
The  pseudorange  that  the  user  measures  is  the  same  regardless  of  which  coordinate  system 
the  user  chooses.  It  is  simpler  to  use  an  ECIF  for  computations  than  to  try  to  make  all 
corrections  in  the  rotating  ECEF  coordinate  system. 

The  fundamental  measurement  of  a  GPS  receiver  is  the  pseudorange.  This  is  a  mea- 
surement of  the  difference  between  the  time  associated  with  the  phase  of  a  GPS  signal 
and  the  time  on  the  local  clock  at  the  arrival  of  that  phase.  We  will  discuss  why  this 
measurement  has  the  same  value  regardless  of  which  coordinate  system  the  user  chooses 
to  do  calculations. 

The  phase  of  an  electromagnetic  wave  is  a  relativistic  scalar.  This  follows  simply 
because  the  phase  is  the  relativistic  scalar  product  of  two  fourvectors:  the  wave  propagation 
fourvector  and  the  position  fourvector.  This  means  the  phase  has  the  same  numerical  value 
in  all  reference  frames.  It  does  not  matter  whether  this  phase  is  evaluated  in  the  earth- 
centered  inertial  frame,  the  rest  frame  of  the  transmitter,  or  the  rest  frame  of  the  receiver. 

For  GPS  signals,  GPS  time  is  encoded  in  the  transmitted  pseudorandom  noise  (PR.N) 
signal.  The  bits  in  the  PRN  code  are  bi-phase  modulated  on  the  L-band  carrier.  Thus 
when  the  time  arrives  for  new  bit  of  the  code,  the  phase  of  the  carrier  either  reverses  180° 
or  it  does  not,  corresponding  to  whether  the  next  bit  is  one  or  zero.  The  receiver  measures 
the  time  of  arrival  of  a  GPS  signal  by  comparing  it  with  a  locally  generated  replica  of  that 
signal.  The  receiver  correlates  the  received  signal  against  the;  local  replica  and  examines 
the  correlator  output  to  find  the  correlation  maximum,  then  locks  the  local  code  replica 
to  the  transmitted  code.  Maintaining  this  lock  in  general  requires  adjusting  the  frequency 
and  time  offset  of  the  code  replica  in  the  correlator  relative  to  the  local  clock.  Tlu^  riMciviM- 
system  then  determines  the  time  difference  between  the  received  phase  tunc  and  the  lime 
on  the  local  clock.   It  is  as  if  the  wave  front  of  tlu^  GPS  signal  luus  an  unambiguous  time 


mai'ker  on  it.  Thus,  the  user  can  measure  the  phase  of  the  received  signal,  hence  the  time 
as  received,  against  the  time  of  the  receiver's  local  oscillator.  The  difference  between  the 
phase  time  as  received  and  the  time  of  the  receiver's  local  oscillator,  after  multiplication 
by  c  to  give  a  quantity  with  units  of  length,  is  the  pseudorange. 

We  can  see  why  the  pseudorange  is  independent  of  which  coordinate  system  the  re- 
ceiver uses  to  solve  for  position  and/or  time.  This  follows  since  the  phase  is  a  relativistic 
scalar  and  the  measurement  is  made  when  the  receiver  clock  is  co-located  with  that  GPS 
signal  phase.  The  measurement  is  essentially  the  time  difference  between  the  local  clock 
time  at  a  particular  point  on  that  clock's  space-time  path,  and  the  GPS  signal  phase  time 
at  the  same  spacetime  point. 

The  navigation  problem  could  be  completely  solved  using  only  the  ECEF  frame.  It  is 
complicated  in  the  ECEF  frame  because  hght  and  hence  GPS  signals  no  longer  travel  in 
Euclidean  straight  lines.  The  time  for  a  signal  to  travel  from  transmitter  at  R  to  receiver 
at  r  is  not  |r  —  R|/c.  Similarly  we  could  set  up  a  coordinate  system  in  either  the  rest  frame 
of  the  receiver  or  of  the  transmitter.  However,  we  emphasize  that  the  problem  is  greatly 
simplified  by  analyzing  the  problem  in  an  ECIF,  then  transforming  results  appropriately 
to  the  ECEF  frame.  This  is  simplified  by  the  fact  that  the  time  coordinate  in  ECIF's 
as  we  have  defined  them  is  identical  to  that  in  the  ECEF  frame.  Thus,  for  example,  the 
relativity  of  simultaneity  is  not  an  issue  because  time,  and  simultaneity,  are  defined  the 
same  way  in  all  these  coordinate  systems. 

2.5  User  Corrections 

At  any  arbitrarily  chosen  instant  the  ECEF  frame  coincides  with  an  ECI  frame  having 
identical  x,  y,  z  axes,  but  not  rotating.  Removing  the  rotation  of  the  earth  from  the  ECEF 
frame  defines  a  coordinate  system  that  is  close  enough  to  an  inertial  frame  to  serve  for 
estimating  the  path  delay.  As  time  passes,  the  ECEF  system  rotates  while  this  ECI  system 
remains  behind,  so  to  speak.  Any  such  fixed  coordinate  system  may  serve  as  an  ECIF  for 
estimating  geometric  range,  using  the  Euclidean  distance  between  the  coordinates  in  the 
ECIF  of  the  satellite  at  transmission  and  the  coordinates  of  the  receiver  at  reception.  The 
time  for  the  signal  to  travel  this  path,  the  geometric  path  delay  tf).  may  be  estimated  to 
high  accuracy  as  the  geometric  range  divided  by  the  defined  speed  of  light  c.  The  slowing  of 
electromagnetic  signal  speed  in  earth's  gravitational  field  (the  Shapiro  time  delay)  amounts 
to  less  than  200  ps  delay  for  GPS  users  near  earth.   [2] 

The  principle  of  the  constancy  of  the  speed  of  light  in  an  inertial  frame  requires  that 
an  ECIF  be  used  for  geometric  path  delay  if  it  is  calculated  by  dividing  the  Euclidean 
distance  by  c.  Using  such  an  ECIF  greatly  simplifies  the  problem  of  solving  for  a  GPS 
user's  position  or  time.  Whereas  GPS  is  intended  to  provide  users  with  their  position 
or  time  in  the  ECEF  system,  if  we  use  ECEF  coordinates  for  geometric  range  with  the 
simple  Euclidean  distance  metric  and  divide  by  c  to  obtain  geometric  path  delay,  we  would 


make  significant  eiiois  since  we  would  be  ignoring  motion  of  the  receiver  during  the  signal 
propagation  time. 

In  general,  a  navigation  solution  requires  pseudorange  measurements  from  at  least 
four  satellites.  These  are  used  to  obtain  geometric  range  delay  estimates,  which  in  turn 
are  used  to  solve  for  position.  To  use  the  simplicity  of  the  Euclidean  distance  metric  that 
comes  with  an  ECIF,  all  satellite  positions  at  the  transmission  epochs  must  be  transformed 
into  the  common  ECIF.  The  user  may  then  solve  for  the  receiver  position  in  this  coordinate 
frame,  and  for  the  GPS  time  t  corresponding  to  this  solution.  Finally,  the  user  must  find 
the  receiver  coordinates  in  the  ECEF  frame  at  the  GPS  time  t  by  accounting  for  the 
rotation  of  the  ECEF  frame  between  the  chosen  moment  tc  at  which  an  ECI  frame  is 
defined,  and  the  GPS  time  t. 

With  a  receiver  using  four  satellites  we  generally  have  five  different  times  of  interest 
as  natural  candidates  for  tc'-  either  a  single  GPS  transmission  time  and  four  different 
reception  times  (time  tagging  at  transmission  time),  or  a  single  reception  time  and  four 
diflPerent  transmission  times  (time  tagging  at  reception  time).  In  any  case  an  ECI  frame 
can  always  be  found  which  coincides  with  the  ECEF  frame  at  a  chosen  instant  of  time,  and 
in  this  sense  the  ECEF  frame  can  determine  an  inertial  coordinate  system.  This  answers 
the  need  for  a  coordinate  system  in  which  light  travels  in  straight  lines  with  speed  c  during 
the  transmission  of  the  signals,  and  in  which  all  the  GPS  satellite  clocks  can  in  principle  be 
synchronized.  We  will  give  a  number  of  examples  below  for  using  this  general  prescription. 
Our  focus  here  will  remain  on  relativistic  corrections  for  GPS  users.  We  do  not  consider 
corrections  for  nonrelativistic  delays:  ionospheric  plasma  delays,  tropospheric  delays  due 
to  water  vapor,  multipath  interference,  or  receiver  system  delays.  Nonrelativistic  delays 
must  certainly  be  accounted  for;  however,  here  our  topic  is  relativity. 

2.6  Relativistic  Doppler  Effect  and  Use  of  the  GPS  Carrier 

Many  receivers  that  make  use  of  the  GPS  carrier  use  only  the  accumulated  phase  to 
smooth  code  measurements.  In  this  case  the  instantaneous  relativistic  Doppler  shift  need 
not  be  considered.  As  we  show  in  Appendix  D,  the  integrated  carrier  as  received  over 
one  cycle  of  the  code  is  equivalent  to  the  change  in  geometric  range  over  that  interval. 
Thus  if  a  receiver  using  accumulated  Doppler  shifts  keeps  track  of  geometric  range,  then 
it  automatically  keeps  track  of  the  relativistic  Doppler  shift.  Nevertheless  some  receivers 
may  use  the  instantaneous  frequency  of  the  received  signals  for  some  purpose.  Such  a 
receiver  must  account  for  the  relativistic  Doppler  effect. 

A  receiver  system  that  uses  the  instantaneous  Doppler  shift  of  the  received  cai'rier 
signal  must  also  correct  for  the  frequency  shift  of  the  received  signal  within  the  framework 
of  relativity.  As  the  true  range  between  the  satelUtc  and  receiver  changes  duo  to  relative 
motion,  the  carrier  frequency  changes  due  to  the  Doppler  effect.  The  relationship  between 
the  received  frequency  /,  and  the  transmitted  frequency,  or  proper  frequency  F.  is  derived 


in  Appendix  D  with  accuracy  sufficient  for  GPS.  This  relationship  is  [2,3,4] 
/  (1  +  A^{rR)/c^  -  v2/2c2)        F{1  +  A^{rT)/c^  -  y^/2c^) 


■.         Nv  ~  1         NV  '  {■^-^^J 


c 


where 

V  =  transmitter  velocity  in  ECI  coordinates, 

V  =  receiver  velocity  in  ECI  coordinates, 

N  =  ECI  unit  vector  in  the  direction  of  signal  propagation, 

A$  =  gravitational  potential  difference  between 

a  clock  and  a  reference  clock  on  the  geoid. 

If  gravitational  potential  differences  can  be  ignored,  the  relationship  between  received  and 
transmitted  frequency  can  be  written  without  approximation  as 

/  F 


7(v)(l-^)        7(V)(1-^)' 
where 


(1) 


7(^)  =      , (2) 


1-    ^ 

The  change  in  frequency  is  closely  related  to  the  rate  of  change  of  range.  In  fact,  in 
Appendix  D  we  show  that  the  relativistic  Doppler  frequency  shift  equation  can  be  derived 
by  differentiating  the  geometric  path  delay  [4,6].  In  Appendix  D  we  have  included  for 
completeness  the  effect  of  the  gravitational  frequency  shift  to  order  c~^,  which  is  sufficient 
for  users  of  the  GPS.  This  proves  the  conceptual  equivalence,  within  the  framework  of 
Special  and  General  Relativity,  of  the  methods  of  pseudorange  (to  obtain  an  instantaneous 
position  solution)  and  integrated  Doppler  frequency  (to  obtain  changes  in  position).  A 
user  of  integrated  Doppler  frequency  accounts  for  changes  in  the  geometric  path  delay 
simply  by  integrating  the  relativistic  Doppler  shift  and  relating  this  observable  quantity 
to  the  change  in  geometric  range.  This  is  discussed  in  detail  in  Appendix  D. 

2.7  The  Rotation  Matrix 

Ahnost  always  the  user  of  the  broadcast  ephemerides  will  want  to  use  coordinate 
transformations  which  correspond  to  rotations  of  the  coordinate  system  about  the  WGS- 
84  z  axis.  For  convenience  we  write  out  here  an  example  of  such  a  rotation  matrix. 
Consider  an  ECIF  with  z  axis  which  coincides  with  the  defined  WGS-84  axis.  Let  the 
position  coordinates  of  some  point  of  interest  in  these  ECI  coordinates  be  denoted  by 


^ECI  = 


xeci 
Veci 
zeci 


(3) 


Suppose  that  the  ECI  axes  coincide  with  the  ECEF  axes  at  the  time  tc-  To  calculate  the 
ECEF  coordinates  of  the  point  at  time  t,  note  that  the  angle  through  which  the  ECEF 
coordinates  have  rotated  is  G{t)  =  fte{t-tc),  where  l^e  =  7.292  115  146  7  x  10~^  rad/s 
is  the  WGS-84  value  of  the  earth  rotation  rate.  Then  the  coordinates  of  the  point  in  the 
ECEF  frame  are  given  by 


XeCEF  =  R    {t)^ECI  — 


cos0(t)       sin0(t)     0 

-sin0(t)     cos0(i)     0 

0  0  1 


XeCI 

Veci 

ZECI 


(4) 


This  is  sometimes  called  a  "passive"  rotation — that  is,  a  rotation  of  the  coordinate  axes 
keeping  physical  position  vectors  unchanged.  In  eq  (4)  the  time  interval  t  —  tc  must  be 
small  since  the  earth  rotation  rate  actually  varies.  With  t  —  tc  as  large  as  3  s,  errors  no 
larger  than  0.05  mm  are  introduced. 

3.  TIME  TAGGING  AT  THE  RECEIVER 

Generally,  a  GPS  navigation  user  measures  the  arrival  times  of  signals  from  four 
different  satellites  and  uses  them  to  solve  for  four  unknowns:  position  x.  y,  z  and  time 
t.  A  common  way  of  doing  so  is  to  make  the  four  measurements  at  one  instant  at  the 
receiver  ("time  tagging  at  the  receiver").  Another  method  is  to  use  signals  which  left 
the  satellites  at  a  common  GPS  time  ("time  tagging  at  the  transmitters").  The  latter 
method  is  more  complex,  both  because  the  signals  are  not  received  simultaneously,  and 
because  GPS  times  on  the  satellite  clocks  must  be  obtained  by  applying  transmitted  clock 
correction  coefficients.  That  is,  the  clocks  on  the  satellites  are  synchronized  only  if  the 
user  applies  transmitted  clock  corrections.  The  GPS  satellite  clocks  themselves  may  differ 
from  GPS  system  time  by  up  to  1  ms  in  the  underlying  ECI  reference  frame  [1.  Sections 
20.3.3.3.1.9-20.3.3.3.3.1].  Thus  in  the  case  of  time  tagging  at  the  transmitters  we  must 
account  both  for  the  motion  of  the  user  during  the  intervals  between  transmission  and 
reception,  and  the  differences  between  SV  clock  time  and  GPS  system  time.  On  the  other 
hand,  as  we  shall  show,  time  tagging  at  the  transmitter  can  be  used  for  navigation  while 
avoiding  most  rotations. 

3.1  Time  Tagging  at  the  Receiver — General  Prescription 

(1)  At  a  chosen  reception  time,  measure  the  transmission  times  using  the  received 
pseudo-random  noise  (PRN)  codes  from  four  satellites.  This  gives  us  the  times 
tsv  according  to  the  SV  clocks  at  transmission  for  the  signals  received  at  the 
chosen  reception  time. 

We  need  the  GPS  time  of  transmission  for  each  satellite  mid  the 

common  reception  time  according  to  the  (possibly  biiuscd)  local  clock. 

We  obtain  the  time  of  reception  from  our  local  clock.  Tlu^  transmission 


time  is  encoded  in  the  received  signal.  Since  we  are  locked  to  the  code,  we 
can  determine  the  offset  of  the  locally  generated  PRN  sequence  required 
to  maintain  lock. 

(2)  Apply  the  prescription  of  Table  20- IV  of  the  ICD-200,  for  corrections  to  obtain 
GPS  system  time  i,  at  transmission. 

This  includes  the  eccentricity  correction  At^,  described  in  the  ICD- 
200,  section  20.3.3.3.3.1,  and  is  an  effect  due  to  relativity.  At^  is  a 
correction  which  applies  to  the  SV  clock;  it  is  the  same  correction  no 
matter  where  the  receiver  is  or  how  the  receiver  is  moving.  The  value  of 
Atr  will  generally  be  different  for  the  clocks  in  different  satellites.  For 
example,  for  a  clock  in  a  satellite  of  orbital  eccentricity  e,  transmitting 
at  an  instant  when  the  eccentric  anomaly  is  E,  the  correction  to  be 
subtracted  from  the  transmitted  time  is 

Atr  =  Fx  10-^°ev^sin  E  ,  (5) 

where 

F  =  -4.442  807  633  s/m^/^  (6) 

A  is  the  semi-major  axis  of  the  satellite  in  meters,  e  is  the  orbital  eccen- 
tricity, and  E  is  the  eccentric  anomaly.  Values  of  these  quantities  are 
transmitted  in  the  navigation  message. 

(3)  Compute  each  satellite's  position  in  the  ECEF  frame  at  its  transmission  time. 

This  is  a  straightforward  application  of  the  equations  of  Table  20-IV 
of  the  ICD-200  [1]  for  determining  the  a:,  y,  z  coordinates  of  the  satellite 
at  the  instant  of  transmission.  Note  that  the  broadcast  message  gives  the 
satellite  positions  at  the  instant  of  transmission  in  ECEF  coordinates — 
specifically  in  the  WGS-84  reference  frame. 

(4)  Choose  an  ECIF  for  computation  of  the  path  delays. 

This  choice  is  arbitrary,  but  some  choices  are  more  convenient  than 
others.  Simplification  may  sometimes  occur  if  the  choice  is  appropriately 
made.  A  natural  choice  for  the  case  of  time  tagging  at  the  receiver  is 
to  freeze  the  ECEF  frame  at  a  time  close  to  the  expected  GPS  time  of 
reception. 

(5)  Transform  the  ECEF  coordinates  of  each  SV  obtained  in  step  (3)  into  the 
chosen  ECIF. 

This  will  normally  require  at  least  three,  and  perhaps  four,  of  such 
ECEF  position  vectors  to  be  rotated.  If  we  freeze  the  ECEF  frame  at 
tlie  instant  of  trHnsmission  of  on(^  of  the  satellites,  only  three  of  these 
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rotations  will  be  needed.  Alternatively,  we  may  freeze  the  ECEF  frame 
at  any  instant  to  define  an  earth-centered  inertial  system.  If  we  choose 
the  GPS  time  equal  to  the  local  clock's  reception  time,  the  later  correc- 
tions may  be  small.  If  the  local  clock  is  sufficiently  close  to  GPS  time 
at  the  instant  of  reception,  then  the  last  step  below,  step  (7),  would  not 
need  to  be  done. 

(6)  Solve  the  path  delay  equations  for  the  receiver's  position  and  time. 

This  can  be  done  by  linearizing  the  propagation  delay  equations 
and  solving  them  iteratively.  Here  we  solve  simultaneously  for  geomet- 
ric range  delay  t]j  and  user  position.  The  user's  motion  through  the 
ECIF  during  propagation  is  included  in  these  to  estimates.  Note  that 
other  contributions  to  path  delay  which  we  do  not  consider  here,  such 
as  ionospheric  and  tropospheric  delays,  should  be  incorporated  in  this 
process.  We  use  an  initial  estimate  of  position  to  solve  for  the  linearized 
corrections  to  receiver  position  and  the  time  offset.  If  these  corrections 
are  not  small  enough  for  the  user,  they  may  be  used  to  obtain  a  next 
estimate  of  position  and  time.  This,  in  turn,  may  be  used  to  obtain 
the  next  linearized  correction  for  position  and  time  offset.  This  iterative 
process  converges  rapidly. 

(7)  Rotate  the  user's  position  coordinates  into  the  ECEF  reference  frame. 

After  finding  the  user's  position  in  the  chosen  ECI  coordinate  sys- 
tem and  the  correct  GPS  system  time,  the  receiver's  ECI  position  co- 
ordinates are  rotated  into  the  ECEF  reference  frame  at  the  instant  cor- 
responding to  the  measured  reception  time.  If  the  ECI  frame  chosen 
in  step  (4)  is  sufficiently  close  to  the  ECEF  frame  at  the  instant  of 
reception,  this  step  may  be  omitted. 

3.2  EXAMPLES 

In  these  examples  we  use  orbital  data  representing  the  GPS  constellation  as  it  ap- 
proximately was  in  March  of  1995.  Values  of  the  semi-major  axes,  eccentricities,  and 
inchnations  have  been  randomized  within  very  small,  but  reasonable  limits.  Keplerian 
orbit  algorithms  are  used  to  propagate  the  SV  orbits  forward  in  time  and  to  accurately 
compute  SV  eccentric  anomalies.  Values  of  Atsvj  are  randomized  using  a  uniform  distri- 
bution within  a  1  ms  range.  The  SV  identification  numbers  used  are  fictitious.  However, 
the  actual  transmitter  and  receiver  data  used  in  these  examples  are  selected  from  among 
many  possibihties  resulting  from  a  simulation  of  the  entire  GPS  constellation.  The  data  are 
thus  obtained  from  a  simulation  where  ''truth"  is  accurately  known  and  can  \yc  c  t)mpai'ed 
with  navigation  solutions. 
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Example  3.21    Mult i- Channel  Receiver,  Time  Tagging  at  the  Receiver 

We  choose  four  satellites  from  this  constellation  in  view  from  our  receiver  location  with 
elevations  above  the  horizon  greater  than  20°.  Suppose  the  receiver  is  truly  at  geocentric 
ECEF  latitude  and  longitude  35°  N,  0°  E,  and  distance  R  =  Q  378  136.3  m  from  earth's 
center.  Suppose  the  user  wants  to  make  a  measurement  of  position  and  GPS  time  at 
his/her  clock  time  t  =  37  239.999  422  365  6  s  of  the  week,  and  that  this  just  happens  to 
be  GPS  time  t  =  37240  seconds  exactly,  for  purposes  of  illustration. 

The  actual  receiver  ECEF  coordinates  at  this  instant,  which  we  give  here  for  compar- 
ison with  the  navigation  solution,  will  be 


■WGS84 


i?  cos  35° 

0 
E  sin  35° 


5  224  663.389  m 

0 
3  658  348.690  m 


(7) 


Such  accuracy  is  not  justified,  but  we  are  specifying  positions  witiiin  1  nun  so  that  the 
convergence  of  the  algorithms  can  be  checked. 

We  follow  the  step  numbering  from  the  prescription  given  above. 

(1)  Suppose  we  already  have  the  measured  SV  clock  times  tsvi  for  the  signals 
received  at  the  chosen  reception  time. 

The  reception  time  is  specified  on  the  local  clock,  which  is  t  = 
37  239.999  422  365  6  s.  This  just  happens  to  be  at  GPS  coordinate  time 
i  =  37  240.000  000  000  0  s  exactly  and  corresponds  to  a  pseudorange 
c  X  0.075  ms  for  SV  1.  The  GPS  time  should  be  recovered  during  the 
computations. 

(2)  We  apply  the  prescription  according  to  Sect.  20.3.3.3.3.1  from  the  ICD- 
GPS-200  describing  the  user  algorithm.  We  then  obtain  the  system  time  t^, 
the  GPS  time  of  transmission,  for  each  satellite  i.  The  relativistic  eccentricity 
correction  is  part  of  this  calculation.  (The  subscript  i  is  not  used  in  the  ICD, 
but  is  added  here  for  clarity.)  The  data  we  actually  use  in  the  present  example 
were  generated  from  a  computer  simulation  which  propagated  orbits  and  clocks 
for  the  entire  GPS  constellation  forward  in  time.  In  this  simulation  the  offsets 
^isvi  were  randomly  generated. 

(3)  Using  the  broadcast  ephemeris,  we  obtain  the  ECEF  coordinates  of  the 
satellites  from  the  prescriptions  given  in  Table  20-IV  of  the  ICD-GPS-200.  From 
the  GPS  time  ti,  the  time  interval  tk  =  ti  —  toe  from  ephemeris  reference  epoch 
is  calculated.  Then  tk  may  then  be  used  in  the  algorithm  for  computation  of 
ephemerides.  This  gives  the  x,  y,  z  coordinates  of  the  sateUite,  in  the  ECEF  frame, 
at  its  transmission  epoch.  Table  1  gives  the  results  after  these  steps.  (These  data 
are  actually  generated  from  a  computer  simulation  of  the  GPS  constellation.) 
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Table  1.  GPS  satellite  positions  in  ECEF  coordinates. 


SV#     Transmission  epoch  ti  {s)  X^  (m)  y^  (m)  Zj  (m) 

1  37  239.924  422  365  6  13  005  878.255  18  996  947.213  13  246  718.721 

2  37  239.920  713  391  8  20  451  225.952  16  359  086.310  -4  436  309.875 

3  37  239.925  307  870  0  20  983  704.633  15  906  974.416  3  486  595.546 

4  37  239.929  346  353  9  13  798  849.321  -8  706  113.822  20  959  777.407 

In  Table  1,  GPS  system  times  ti  rather  than  the  time  intervals 
tk  —  ti  —  toe  from  ephemeris  reference  epoch,  are  used  to  label  the 
events.  The  subscript  i  varying  from  1  to  4  labels  data  from  the  different 
satellites.  Table  1  assumes  that  the  eccentricity  correction  discussed  in 
Appendix  A  has  already  been  applied  to  obtain  the  transmission  epoch 
data  given  in  the  second  column  of  the  table.  The  position  coordinates 
of  each  satellite  are  given  in  the  ECEF  frame,  which  is  turning  in  inertial 
space. 

(4)  To  rotate  from  ECEF  to  ECIF  coordinates,  the  inverse  of  the  rotation 
matrix  of  eq  (4)  is  required.  Thus  we  need  {R^{t^  —  ic))~\  which  must  be 
calculated  and  applied  to  the  coordinates  of  each  satellite  individually.  That  is, 
ti  will  be  different  for  each  satellite,  but  tc  is  the  same  for  all  satellites.  For 
this  example,  we  choose  an  ECIF  by  fixing  the  rotation  of  the  ECEF  frame  at 
tc  =  37  239.000  000  000  0  s.  We  stress  that  this  choice  is  arbitrary — better  choices 
are  usually  available.  The  purpose  of  this  choice  is  to  illustrate  its  arbitrariness 
and  also  to  illustrate  how  rapidly  the  iteration  algorithm  converges.  With  this 
choice,  for  example,  the  rotation  matrix  for  satellite  3  is 


{R^{t3-tc)y 


cosCte{t3  -tc)     -sinpeitd  -tc)     0 

sintteitz  -tc)       cos Clf,{t3  -tc)       0 

0  0  1 


and  fteits  -  tc)  =  0.925  307  870  0  ■  fte  =  6.747  451  534  x  10"^  Corresponding 
transformations  are  applied  to  each  of  the  transmitted  data  sets. 

(5)  We  transform  the  coordinates  in  Table  1  into  this  ECIF.  The  rotation 
matrix  is  different  for  each  transmission  epoch.  We  use  a  subscript  C  to  indicate 
that  the  inertial  system  has  been  arbitrarily  chosen.  The  results  from  transform- 
ing to  the  ECI  frame  are  given  in  Table  2. 

(6)  The  propagation  delay  equations  to  be  solved  arc  discussed  in  Api)cndi,\ 
B.  We  give  a  quick  summary  here.  The  notation  is  cis  follows.  The  four  GPS 
satellites,  at  GPS  times  ti,  t2,  ^3,  t^  send  out  signals  from  the  ECIF  locations 
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(8) 


Table  2.  Transmitted  data  transformed  to  the  chosen  ECI  system 

SV  #  Transmission  Epoch  t^  (s)  Xa  (m)  ya  (m)  ZCi  (ni) 

1  37  239.924  422  365  6  13  004  597.642  18  997  823.895  13  246  718.721 

2  37  239.920  713  391  8  20  450  127.566  16  360  459.358  -4  436  309.875 

3  37  239.925  307  870  0  20  982  631.270  15  908  390.245  3  486  595.546 

4  37  239.929  346  353  9  13  799  439.294  -8  705  178.668  20  959  777.407 


ri5  1*2  5  1*3  5  r4  5  respectively.  (The  subscripts  now  label  data  from  the  different 
satellites,  after  transforming  to  the  chosen  ECI  reference  frame.  The  rotations  do 
not  affect  the  times.)  These  four  signals  are  received  simultaneously  at  time  i  by  a 
receiver  at  position  re-  The  problem  is  to  determine  t  and  r^  at  the  receiver.  The 
velocity  of  the  receiver  does  not  have  to  be  considered.  The  receiver's  position 
re  at  the  time  of  reception  t  will  be  determined  by  the  solution  of  four  equations 
which  express  the  condition  that  the  speed  of  propagation  is  c.  The  four  equations 
to  be  solved  are 

(re  -  Tjf  -  c\t  -  tjf  =0;     3  =  1,2,3,4.  (9) 

However,  if  the  position  yc  and  time  t  are  known  approximately,  the  equations  can 
be  reduced  to  a  system  of  four  linear  equations  which  can  be  solved  by  standard 
matrix  inversion  techniques.  These  equations  are  derived  in  Appendix  B,  and  are 

(r(^)  _  r^.)  .  Ar  -  c^{t^'^  -  tj)At  =  -  (^c^{t'^'^  -  tjf  -  (r^''^  -  r,)^)  .  (10) 

In  eq  (10),  the  quantity  c{t^^^  —tj)  is  the  ith  estimate  of  the  pseudorange  from  the 

(i) 

receiver  to  the  jth  satellite,  and  r^  is  the  zth  estimate  of  the  receiver  position  in 
ECI  coordinates.  Equation  (10)  is  a  system  of  Hnear  inhomogeneous  equations 
in  the  corrections  Ar,  At,  which  can  be  solved  by  matrix  inversion.  The  matrix 
of  coefficients  of  the  unknowns  Ar  and  At  will  usually  be  nonsingular,  unless  the 
configuration  of  satellites  is  so  unfavorable  that  eqs  (9)  do  not  have  a  solution 
(such  a  condition  is  possible.)  The  solutions  obtained  will  be  approximate,  but 
can  be  used  to  obtain  new  trial  values;  the  iteration  can  be  repeated  as  many 
times  as  necessary  to  obtain  the  accuracy  required. 

To  illustrate  this  process  in  the  present  case,  as  our  initial  guess  at 
receiver  position  we  take  a  worst  case  and  assume  the  receiver  is  at  the 
center  of  the  earth.  Also,  we  do  not  know  the  receiver  clock  bias.  In 
this  example  we  assume  the  clock  bias  is  such  that  the  reception  event 
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occurs  at  t  =  ti  +0.07500  s.  Then  Table  3  gives  the  results  at  each  stage 
of  the  iteration. 

The  calculation  converges  to  within  1  m  after  three  iterations.  The 
time  at  reception  on  the  receiver  clock,  determined  by  the  solution,  is 

tR  =  37  240.000  000  000  Os . 


(7)  These  results  must  lastly  be  transformed  into  the  WGS84  system  by 
applying  the  rotation  R^itu  —  tc)-  That  is,  to  obtain  WGS84  coordinates  at 
the  instant  of  reception,  we  must  use  a  sidereal  rotation  corresponding  to  that 
instant.  The  necessary  rotation  matrix  is 


R^tR-tc)  = 


COSfteitR  -tc) 

-  sin  ft,,{tR  -  tc] 
0 


sniVt^{tR-tc)      0 

COS^e{tR-tc)       0 

0  1 


(11) 


and  Cte{tR  —  tc)  —  7.292  115  146  7  x  10~^  because  in  this  case,  tR  —  tc  happens 
to  be  exactly  1  s.  Upon  applying  this  rotation  to  the  position  coordinates  given 
in  the  last  line  of  Table  3,  the  measured  values  of  the  receiver  position  are 


■WGSBA 


5  224  663.388  m 

•  0.000  m 
3  658  348.689  m 


(12) 


These  agree  with  the  true  position  coordinates  to  within  a  inillimeter. 


Table  3.  Results  of  iterative  solution  of  propagation  delay  equations. 


Trial 

# 

User's  X  position 

User's  y  position 

User's  z  position 

User's  clock  time 

0  (start) 

0  m 

0  m 

0  m 

37  239.999  422  365  6  s 

1 

5  057  363.392  m 

2  355.126  m 

3  541  092.792  m 

37  239.997  532  305  8  s 

2 

5  226  931.552  m 

354.224  m 

3  659  938.391  m 

37  240.000  033  455  9  s 

3 

5  224  663.780  m 

380.983  m 

3  658  348.973  ui 

37  240.000  000  006  0  s 

4 

5  224  663.374  m 

380.988  m 

3  658  348.689  ni 

37  240.000  000  000  0  s 

5 

5  224  663.374  m 

380.988  m 

3  658  348.689  m 

37  240.000  000  000  (1  s 
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Table  4.  Iterative  solution  ignoring  earth  rotation. 


0  1 

(staxt) 

0 

0 

0 

1 

5  057  359.485 

2  012.579 

3  541  090.386  m 

2 

5  226  927.086 

00.194 

3  659  935.704  m 

3 

5  224  659.325 

27.107 

3  658  346.291  m 

4 

5  224  658.919 

27.112 

3  658  346.008  m 

5 

5  224  658.919 

27.112 

3  658  346.008  m 

Trial  #        User's  x  position  (m)     User's  y  position  (ni)     User's  z  position  (m)     User's  clock  time  (s) 

37  239.999  422  365  6 
37  239.997  532  323  0 
37  240.000  033  469  3 
37  240.000  000  019  6 
37  240.000  000  013  6 
37  240.000  000  013  6 


Example  3.211    Erroneous  Use  of  ECEF  Coordinates 

One  might  think  that  effects  due  to  earth  rotation,  during  the  very  small  signal  prop- 
agation times  from  transmitters  to  receiver,  were  sufficiently  small  that  transforming  to 
ECI  coordinates  was  not  necessary.  Let  us  therefore  test  this  with  the  example  data  given 
in  Table  1.  We  again  initialize  the  iteration  as  in  Example  3.21,  above.  The  results  of  the 
iteration,  using  the  same  algorithm  as  before  are  given  in  Table  4. 

In  this  case,  the  solutions  are  in  error  by  a  significant  amount.  The  receiver  time  is 
off  by  14  ns  and  the  position  is  in  error  by  almost  30  m.  Thus,  it  is  an  error  to  ignore 
earth  rotation  during  the  propagation  of  signals  from  transmitters  to  receiver. 

However,  these  results  suggest  an  approach  which  could  result  in  a  computationally  ef- 
ficient strategy  for  the  receiver.  Suppose  that  the  receiver  is  continually  monitoring  signals 
from  four  or  more  satellites  and  keeping  track  of  its  approximate  position  by  neglecting 
earth  rotation,  as  in  Example  3.211  above,  but  that  an  accurate  position  is  not  needed 
until  some  particular  moment.  At  that  moment,  the  erroneous  position  and  time  could  be 
used  as  a  starting  point  for  a  correct  iterative  solution  of  the  one-way  light  time  equations. 
Since  the  starting  point  is  not  far  off,  this  iteration  should  converge  very  rapidly.  Thus, 
transformations  to  the  chosen  inertial  frame  C  would  be  made  only  when  an  accurate 
position  and  time  were  needed,  eliminating  much  computation  at  other  times. 

Example  3.212    Iteration  with  a  Succession  of  Inertial  Frames 

Another  approach  is  to  consider  a  succession  of  inertial  frames  in  the  iteration  process, 
each  of  which  is  aligned  with  the  WGS84  frame  at  the  assumed  reception  time.  This  would 
have  the  nice  advantage  that  at  the  end  of  the  iteration,  the  ECIF  would  be  exactly  ahgned 
with  the  WGS84  frame,  so  no  final  rotation  would  be  necessary.  The  disadvantage  is  that 
the  ECIF  coordinates  of  the  transmission  events  must  be  recomputed  at  each  stage  of 
the  iteration.  Thus  in  the  intermediate  stages  of  iteration  more  computation  is  required. 
During  actual  use,  we  might  have  very  good  estimates  of  the  receiver  position  available 
to  start  the  iteration.  Thus,  very  few  iterations  could  yield  receiver  position  to  sufficient 
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accuracy,  making  this  approach  favorable. 

To  illustrate  this  approach,  we  shall  use  a  succession  of  equations  of  the  form 

Xc  =  R^{tc){R^r'^WGSS4  -  {R^{t  -  tc)r'XwGS84  (13) 

to  obtain  inertial  coordinates  of  the  transmission  events.  We  will  first  make  a  reasoned 
guess  for  the  receiver  time  tn  and  then  set  tc  —  iR-  One  iteration  will  yield  a  more  accurate 
estimate  for  Ir.  We  will  then  choose  tc  to  be  this  new  value  and  recompute  satellite 
transmission  event  data  for  this  new  tc-,  and  repeat  the  process  as  often  as  necessary.  The 
results  are  given  in  Table  5. 

Example  3.213    The  Benefit  of  Better  Initialization 

It  pays  to  be  more  careful  with  our  initial  guesses.  To  illustrate  this,  let  us  take  as 
our  initial  time  the  latest  transmission  epoch  plus  75  ms  and  take  the  initial  position  to  be 
the  subsatellite  point  of  the  latest  transmission  epoch.  Again  using  a  succession  of  inertial 
frames  as  in  the  above  example,  the  results  are  given  in  Table  6.  What  has  been  achieved 
is  more  rapid  convergence,  compared  with  the  previous  example. 

Table  5.  Iterative  solution  with  a  succession  of  inertial  frames. 

Trial  #       User's  x  position  (m)   User's  y  position  (m)     User's  z  position  (m)   User's  clock  time  tc  (s) 

37  239.999  422  365  6 
37  239.997  532  305  8 
37  240.000  033  455  9 
37  240.000  000  006  0 
37  240.000  000  000  0 


Table  6.  Effect  of  better  initialization. 


0  (start) 

0 

0 

0 

1 

5  057  363.550 

1  986.551 

3  541  092.792 

2 

5  226  931.564 

-25.990 

3  659  938.391 

3 

5  224  663.794 

-0.018 

3  658  348.973 

4 

5  224  663.388 

0.000 

3  658  348.689 

Trial 

# 

Ui 

ser's  X  position  (m) 

Ui 

ser's  y  position  (m) 

Ui 

ser's  z  position  (m) 

User's  clock  time  tc  (s) 

0  (start) 

3  313  469.280 

-2  090  568.570 

5  032  997.819 

37  240.004  346  353  9 

1 

5  190  759.859 

400.888 

3  634  586.629 

37  239.999  499  918  6 

2 

5  224  663.389 

-0.892 

3  658  412.593 

37  240.000  000  000  0 

3 

5  224  663.389 

-0.000 

3  658  348.689 

37  240.000  000  000  0 

4 

5  224  663.389 

0.000 

3  658  348.689 

37  240.000  000  000  0 
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4.  TIME  TAGGING  AT  THE  TRANSMITTER 

Let  us  now  consider  time  tagging  at  the  transmitters  with  four  sateUites.  This  ap- 
proach requires  a  prior  estimate  of  receiver  velocity  and  perhaps  acceleration  over  the 
interval  of  signal  reception.  If  all  signals  leave  satellites  simultaneously  with  respect  to  the 
ECEF  frame,  most  users  will  receive  these  signals  within  20  ms  of  each  other.  For  most 
users  acceleration  will  contribute  negligibly. 

We  address  two  separate  issues  for  this  type  of  receiver.  The  most  important  from 
a  conceptual  point  of  view  is  that,  during  the  20  ms  or  so  which  may  be  required  for 
the  receiver  to  make  four  time-of-arrival  measurements,  the  receiver  will  generally  move 
through  the  inertial  frame.  This  motion  must  be  accounted  for  in  order  to  achieve  accurate 
position  determination. 

A  second,  more  practical  comphcation,  is  that  the  C/A  code  repeats  every  millisecond, 
and  it  might  be  easiest  to  choose  the  GPS  system  time  tag  tc  for  the  measurements  to 
be  identical  to  one  of  the  instants  at  which  the  C/A  code  begins  to  repeat.  We  call  this 
a  "code  boundary."  In  the  present  example  we  shall  illustrate  use  of  C/A  code  in  which 
the  receiver  only  measures  the  arrival  times  of  code  boundaries  on  its  local  clock.  We 
shall  first  discuss  C/A  code  boundary  measurements  and  the  necessary  time  measurement 
corrections.  Then  we  shall  discuss  how  to  account  for  receiver  motion. 

4.1  Code  Boundary  Measurements 

First  we  discuss  the  choice  of  the  time  tc  at  which  the  ECEF  frame  is  to  be  "frozen" 
in  place  to  create  the  ECIF.  The  time  tc  denotes  a  GPS  system  time.  Unfortunately, 
because  of  the  clock  and  relativistic  corrections  represented  by  the  quantity  Atsv  in  the 
ICD-GPS-200,  the  GPS  system  time  t  at  a  satellite  differs  from  the  satellite  clock  time, 
tsv  according  to  the  equation 

t^tsv-Atsv-  (14) 

Thus,  transmission  of  a  a  particular  system  time  tc  would  occur  at  satellite  clock  time 
which  is  numerically  Atsv  later  (or  earlier  if  Atsv  is  negative)  by 

tsv  ^tc  + Atsv.  (15) 

The  correction  Atsvj  foi'  a  particular  satellite  is  contained  in  the  navigation  message 
and  can  be  computed  by  the  receiver.  It  changes  relatively  slowly  during  the  20  ms  or  less 
required  by  the  receiver  to  make  four  time-of-arrival  measurements,  so  we  shall  assume 
that  for  a  particular  set  of  four  time-of-arrival  measurements,  Atsvj  are  constants. 

To  assist  in  making  clear  the  distinctions  among  the  various  times  which  occur  in  the 
following  discussion,  we  shall  use  a  subscript  R  for  times  measured  on  the  receiver  clock, 
and  SV  for  times  measured  on  the  transmitter  clock.  Both  of  these  may  differ  from  GPS 
system  times  for  events,  due  to  noise,  drifts,  offsets,  and  relativistic  effects. 
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tp  on  SV  clock 


GPS  system  time 

tQ  =  chosen  GPS  system  time 
Figure  1.  Illustration  of  SV  and  receiver  clock  offsets. 

4.2  Phase  Time 

The  "phase  time"  r  is  transmitted  from  the  satellite  at  SV  clock  time  tsv  and  propa- 
gates to  the  receiver  where  it  arrives  at  receiver  clock  time  Ir.  In  a  particular  measurement, 
r  is  some  number.  If  this  number  is  an  integral  multiple  of  1  ms  so  that  it  lies  on  a  code 
boundary,  then  it  may  be  particularly  convenient  for  the  receiver  to  implement  the  mea- 
surement of  the  arrival  time.  We  shall  assume  in  the  following  that  r  lies  on  a  code 
boundary.  Thus  for  the  channel  corresponding  to  the  jth  satelhte  (with  j  =  1,  2,  3,  4),  the 
sequence  of  events  may  be  described  as  follows. 

(a)  A  common  system  time  tc  is  chosen  for  the  measurements,  which  are  to 
be  time  tagged  at  the  transmitter  at  this  system  time  tc-  This  system  time  tc 
is  some  integral  multiple  of  1  ms  and  thus  there  is  a  code  boundary  in  the  phase 
time  which  is  numerically  equal,  t  =  tc- 

(b)  While  the  phase  code  is  numerically  equal  to  tc  and  is  transmitted  when 
the  SV  clock  time  tsv 3  =  t  =  tc-,  this  does  not  usually  occur  at  the  GPS  system 
time  t  =  tc-  Because  of  clock  and  relativistic  corrections,  transmission  of  the 
true  GPS  system  time  tc  occurs  ^tsvj  later. 

(c)  The  wave  train  propagates  to  the  receiver  where  the  code  boundary  r  = 
tc  is  detected  and  the  receiver  clock  time  tjij  of  this  event  is  recorded. 

This  sequence  of  events  is  illustrated  in  Figure  1. 


The  arrival  of  a  point  in  the  code  corresponding  to  tru(>  GPS  system  tinu^  ic  ^vill 
occur  an  amount  A^;?j  later.  Because  of  the  relative  motion  of  the  riHcMVtn  and  tiansniitlcr 
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through  the  ECI,  Atjij  will  in  general  be  different  from  Atsvj  due  to  the  first-order  Doppler 
effect.  The  receiver  needs  to  estimate  this  time  difference  in  order  to  find  the  arrival  time, 
on  the  local  clock,  of  the  point  in  the  code  corresponding  to  the  true  time  tag  tc- 

The  small  interval  Atjij  can  be  very  accurately  estimated  by  using  the  Doppler  for- 
mula, eq  (D.6)  in  Appendix  D,  which  relates  the  coordinate  time  interval  dtn  between 
signal  arrivals  at  the  receiver,  and  the  corresponding  time  interval  dtsv  at  the  transmitter 

/        N-v\                 /        N-V\ 
dtRil =  dtsv     1 ,  {D.6) 


where  V  and  v  are  the  transmitter  and  receiver  velocities  through  the  ECI,  respectively, 
at  the  instants  of  transmission  and  reception.  This  gives  the  relation  between  the  time 
interval  separating  transmission  of  two  nearly  equal  phase  times,  and  the  time  interval 
separating  reception  of  the  two  phase  times.  The  time  intervals  dtsvj  are  less  than  1  ms. 
The  first-order  Doppler  corrections  are  typically  of  order  10~^.  Therefore  we  need  good 
estimates  of  the  velocities  of  receiver  and  transmitter,  so  that  the  Doppler  corrections  are 
known  within  10~^  in  order  to  calculate  dtn  to  better  than  1  ns. 

Because  the  Doppler  corrections  in  eq  (D.6)  are  small,  we  can  use  eq  (D.6)  for  small 
time  intervals  without  making  appreciable  error.  Thus,  to  sufficient  accuracy  for  GPS,  eq 
(D.6)  can  be  rewritten  as 

Atn,=Atsvj^r^-^^.  (16) 

The  receiver  needs  estimates  of  the  receiver  velocity  v  and  the  unit  vector  N  pointing 
from  transmitter  to  receiver  in  the  ECI.  These  can  be  obtained  from  estimates  of  position 
and  velocity  which  have  been  built  up  over  time  during  the  computation  of  successive 
navigation  solutions.  The  transmitter  velocity  Vj  can  be  computed  from  the  prescriptions 
given  in  Table  20-IV  in  the  ICD-GPS-200,  by  taking  derivatives  of  the  satellite  position 
with  respect  to  coordinate  time.  Since  all  these  vectors  must  be  represented  in  the  chosen 
ECI,  some  rotations  may  be  necessary.  It  is  useful  to  keep  in  mind,  however,  that  dot 
products  may  be  computed  in  any  convenient  ECI  frame. 

Thus,  in  brief,  if  the  phase  time  corresponding  to  true  GPS  system  time  tc  is  later 
than  tsvj  —  tc  by  Atsvj,  then  the  time  of  reception  of  the  phase  time  corresponding  to 
tc  at  the  receiver  is  later  by  Atjij  given  by  eq  (16)  above.  These  time  relationships  are 
illustrated  in  Figure  1. 

4.3  Accounting  for  Receiver  Motion 

In  the  following  discussion  we  will  assume  that  the  code  boundary  measurements  have 
been  made  and  corrected  as  described  above.  After  making  the  necessary  corrections, 
denote  the  values  of  the  arrival  times  of  the  code  corresponding  to  the  true  GPS  time 
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tc  on  the  receiver  clock  by  tcj  =  tRj  +  ^tjij.  These  will  still  not  be  GPS  system  times 
because  in  general  there  will  be  a  bias  on  the  local  clock.  Denote  the  bias  by  be,  and 
denote  the  GPS  system  times  of  arrival  by  tj;  then  the  relation 

tj=tcj+bc  (17) 

holds.  Note  that  there  is  no  subscript  R  on  the  time  on  the  left  side  of  eq  (17)  because 
this  symbol  now  refers  to  the  GPS  system  time  of  an  event. 

Further,  because  of  the  receiver's  motion  through  the  ECI,  the  receiver  will  not  be  at 
a  fixed  position  during  the  measurements,  which  occur  from  67  to  86  ms  after  tc  if  the 
receiver  is  near  the  earth.  This  means  it  is  necessary  to  have  available  to  the  receiver  an 
estimate  of  receiver  velocity  in  the  chosen  ECI  frame,  so  that  the  receiver  positions  at  the 
various  reception  times  can  be  calculated.  Let  the  receiver  position  and  velocity  in  the  ECI 
at  the  system  time  tc  be  re  and  vc  respectively.  We  assume  vc  is  known  to  sufficient 
accuracy,  but  re  and  be  are  to  be  determined  by  the  measurement.  Thus  the  positions  of 
the  receiver  at  the  measurement  times  are  given  by 

rRj  =  TRitj)  =  rc  +  ^citcj  +bc-  tc) .  (18) 

The  sateUite  positions  at  the  instant  of  transmission  of  the  signals  at  the  common  GPS 
system  time  tc  wiU  be  denoted  by  r^.  These  positions  are  obtained  from  the  navigation 
messages  and  are  specified  in  the  WGS-84  reference  frame.  This  means  that  all  of  these 
position  vectors  are  specified  in  the  same  reference  frame.  Therefore  it  is  natural  when 
setting  up  the  propagation  delay  equations  to  choose  an  ECI  reference  frame  which  coin- 
cides with  the  ECEF  frame  at  the  instant  tc-  Thus  no  rotations  are  necessary  in  order 
to  express  these  SV  positions  in  a  common  inertial  frame.  If  we  solve  for  re,  this  will  be 
in  the  ECIF  (and  the  ECEF  frame)  at  tc,  so  no  rotation  will  be  required  for  that  vector. 
However,  if  we  want  trj  in  the  ECEF  frame,  rotations  will  be  necessary. 
The  propagation  delay  equations  are 

i^Rj  -^jf  =  c^(^j  -tcf 


or 


(re  +  ^citcj  +bc-  tc)  -  Fjf  -  c^itcj  +bc-  tc)^  .      j  =  1.2.3.4.  (19) 


These  four  equations  can  be  solved  iteratively  for  the  unknown  bias  be  and  the  tlu-ee 
position  components  re  of  the  receiver  in  the  ECEF  frame  at  the  chosen  system  time  tc- 
One  iterative  method  of  solution  is  described  in  Appendix  C.  We  briefly  summai-ize  that 
process  here. 

We  write  the  position  and  bias  in  terms  of  estimates,  plus  small  corrections: 

re  =  r^c  +  ^r-      ^c  =  i^c  +  ^^  ■  (20) 
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Substituting  these  expressions  into  eq  (19)  and  expanding,  keeping  only  linear  terms  in  the 
small  quantities  Ar  and  At,  leads  to  a  linear  system  of  four  equations  in  four  unknowns: 


R^;^  -Ar-c^iT 


(7^^)-v..Ry))A.4    c^(lj^O'-K) 


(21) 


where 


and 


T^/^=tcj+b)^^-tc,  (22) 


Rf  =r^)+vcTf)-r,.  (23) 

This  system  of  equations  is  solved  by  iteration  in  the  manner  described  in  Appendix  C. 

4.4  Time  Tagging  at  the  Transmitter — General  Prescription: 

We  now  describe  the  specific  procedure  in  parallel  with  the  steps  described  above  for 
time  tagging  at  the  receiver. 

(1)  Measure  the  times  tjij  on  the  local  clock  at  the  reception  of  a  given  code 
epoch  r,  corresponding  to  a  single  GPS  phase  time  r  in  each  code  train,  from 
each  of  four  satellites. 

It  may  be  convenient  for  a  receiver  to  measure  the  time  on  its  local 
clock  for  the  arrival  of  a  particular  point  of  the  received  code  train, 
such  as  one  of  the  code  boundaries.  A  code  boundary  corresponds  to  a 
particular  millisecond  of  the  week.  This  millisecond  is  our  chosen  time 
tag.  However,  each  SV  clock,  and  hence  the  clock's  transmitted  code, 
is  generally  offset  from  from  GPS  time,  and  the  offset  will  be  different 
for  each  clock.  The  differences  will  be  less  than  1  ms.  We  will  need  to 
account  for  this  in  step  2. 

(2)  Apply  the  prescription  for  corrections  as  given  in  the  "User  Algorithm  for  SV 
Clock  Correction"  of  the  ICD-200,  section  20.3.3.3.3.1,  to  obtain  Atsvj 

For  a  given  satellite,  Atsv  is  the  difference  between  the  instant  the 
code  for  phase  time  r  was  transmitted  and  the  value  of  GPS  time  r. 
(See  eqs  (14-15).)  Atsv  includes  both  the  relativistic  correction  Atr 
and  the  broadcast  second-order  polynomial  correction  for  the  SV  clock 
offset  from  GPS  time. 

(3)  Compute  each  satellite's  position  in  the  ECEF  frame  at  the  time  of  transmis- 
sion of  GPS  time  t. 

This  is  a  straightforward  application  of  the  equations  of  Table  20- 
IV,  for  determining  the  .x,  t/,  z  coordinates  of  the  satellite  at  the  GPS 
time  tc- 
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(4)  Choose  an  ECIF  frame  for  computation  of  the  one-way  hght  times. 

In  this  case  the  natural  choice  is  to  freeze  the  ECEF  frame  at  the 
GPS  transmission  time  tc- 

(5)  If  the  ECI  is  defined  by  freezing  the  ECEF  frame  at  GPS  transmission  time 
tc,  the  coordinates  of  each  sateUite  will  not  need  to  be  rotated.  If  we  estimate 
receiver  position  at  tc  we  do  not  need  to  rotate  the  receiver  position  coordinates 
either.  We  need  to  adjust  the  local  clock  measurements,  however,  for  the  time 
interval  Atsv  and  the  motion  of  the  receiver  over  the  time  Atsv-  This  is  because 
the  code  boundary  numerically  corresponding  to  tc  is  not  necessarily  transmitted 
at  GPS  time  tc- 

We  may  use  a  deterministic  estimate  of  the  change  in  receiver  po- 
sition over  the  reception  interval.  Since  this  model  is  usually  a  velocity, 
we  will  refer  to  it  that  way  here  for  simplicity.  This  velocity  estimate 
remains  fixed  during  the  solution  of  user  position  and  time.  It  is  neces- 
sary to  obtain  the  velocity  estimate  in  the  coordinates  of  the  chosen  ECI. 
Then  we  may  estimate  receiver  coordinates  in  the  ECEF  frame  which 
coincides  at  tc  with  the  ECI.  We  then  use  the  velocity  to  estimate  re- 
ceiver positions  at  the  other  reception  times.  We  assume  the  user  clock 
offset  from  GPS  time  is  constant  over  the  interval  of  reception. 

First  apply  the  corrections  represented  in  eq  (16)  to  compute  Atfij. 
Then  compute  the  times  of  arrival  on  the  local  clock  tcj,  of  the  signal 
transmitted  at  GPS  time  r.  Iftjij  is  the  local  clock  time  for  the  received 
epoch  from  satellite  j,  then 

tcj  =  tRj  +  AtRj  (24) 

is  the  local  clock  time  of  reception  of  the  signal  transmitted  at  GPS  time 
t  from  satellite  j.  The  correction  AtRj  is  calculated  from  eq  (16)  using 
estimates  of  vrj,  Tj,  and  Vj,  and  Vj  in  the  ECI  which  have  been  built 
up  from  a  filtering  algorithm  or  from  prior  solutions. 

The  time  relationships  are  illustrated  in  Figure  1 .  For  a  given  satel- 
lite vehicle  SV,  Atsv  is  the  diff'erence  between  the  instant  the  code  for 
GPS  time  r  was  transmitted  and  the  value  of  GPS  time  r.  Atsv  includes 
both  the  relativistic  correction  Atr  and  the  broadcast  second-order  poly- 
nomial correction  for  the  SV  clock  offset  from  GPS  time.  Atjij  is  related 
to  Atsv  by  first-order  Doppler  corrections  given  in  eq  (16). 

(6)  Solve  the  equations  for  the  geometric  path  delays  to  find  the  receiver's  position 
and  time  offset  at  tc- 

Again  we  linearize  the  propagation  delay  equations  and  solve  Ihcm 
iteratively.    We  may  iterate  to  find  the  position  at  the  first  reception 
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time,  using  the  velocity  estimate  to  extend  to  the  other  reception  times. 
This  iteration  is  similar  to  step  (6)  in  the  case  of  time  tagging  at  the 
receiver. 

(7)  If  we  have  solved  for  the  user  position  at  tc  there  is  no  final  rotation.  In 
this  case  the  ECIF  and  ECEF  frame  coincide  at  tc-  If  we  want  the  ECEF  frame 
position  at  a  different  reception  time,  a  rotation  will  be  required. 

It  would  be  possible  to  update  the  velocity  estimate  using  positions  obtained  from 
two  transmission  time  tags.  This  leads  into  other  concerns  and  options  associated  with 
using  GPS:  filtering  estimates  over  time.  Techniques  for  estimating  position,  velocity,  and 
acceleration  could  be  coupled  with  strategies  for  filtering  these  estimates  over  time.  We 
will  not  discuss  these  options  here.  We  will  mention,  however,  that  if  there  are  infrequent 
measurements  of  GPS  signals  from  individual  satellites,  the  local  clock  may  be  used  to 
flywheel  between  them  and  find  a  solution.  This  may  require  careful  filtering  algorithms 
for  estimation  of  the  local  clock  frequency.  In  particular,  relativistic  effects  on  the  local 
clock  frequency  due  to  velocity  or  gravitational  potential  may  have  to  be  considered. 

4.5  EXAMPLES 

Example  4.51  Time  Tagging  at  Transmitter 

We  give  here  an  example  of  the  analysis  of  data  when  the  measurements  are  time 
tagged  at  the  transmitter.  This  means  that  the  receiver  measures  the  arrival  times  of  at 
least  four  signals  which  were  transmitted  simultaneously,  at  GPS  system  time  tc,  from 
four  satellites  which  are  labelled  with  index  j  =  1, 2, 3, 4.  This  method  of  time  tagging  has 
the  great  advantage  that  the  satellite  positions  at  the  instant  of  transmission,  as  specified 
in  the  navigation  message,  are  all  given  in  the  same  coordinate  system.  Thus  no  rotations 
are  necessary  to  express  these  positions  in  a  common  inertial  frame.  For  this  example  we 
choose  a  code  boundary  exactly  corresponding  to 

tc  =  240.000  000  000  0  s . 

The  following  steps  are  in  the  order  of  the  seven  steps  listed  in  Sections  2  and  3  above. 

(1)  The  times  tfij  on  the  local  clock  of  reception  of  the  code  epoch  r  =  240 
sec  exactly  are  recorded.  These  are  given  in  the  first  column  of  Table  7. 

(2)  The  satellite  clock  corrections  Atsv,  including  the  relativity  corrections, 
are  obtained  from  the  navigation  messages;  these  are  given  in  the  second  column 
of  Table  7.  In  this  simulation  these  offsets  are  randomly  generated  within  a  1  ms 
range  using  a  uniform  distribution. 

(3)  Compute  each  satellite's  position  in  the  ECEF  frame  at  the  time  of 
transmission  of  GPS  system  time  r. 
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This  is  a  straightforward  appHcation  of  the  equations  of  Table  20-IV 
in  the  ICD-200,  for  determining  the  x,  y,  z  coordinates  of  the  satelhtes  at 
the  GPS  time  tc-  The  results  are  given  in  Table  7  for  the  four  satelUtes. 

(4)  Choose  an  ECI  frame  for  computation  of  the  one-way  Hght  times. 

In  this  case  the  natural  choice  is  to  freeze  the  ECEF  frame  at  the 
GPS  transmission  time  tc- 

(5)  Compute  the  received  time  offsets  At/^j  and  the  times  tcj  of  arrival  on 
the  local  clock,  of  the  signals  transmitted  at  GPS  system  time  tc  =  t. 

We  must  make  the  first-order  Doppler  corrections,  then  compute 
the  time  intervals  ^tnj.  This  requires  estimating  the  Doppler  correc- 
tion. However  since  the  Doppler  corrections  are  very  small,  in  a  first 
approximation  we  might  neglect  the  Doppler  correction  factor  entirely 
and  assume  Ai/?j  ■^  Atsvj-  We  might  then  use  the  resultant  positions 
to  estimate  Nj  to  within  a  few  percent.  If  the  estimate  of  position  in  the 
chosen  ECIF  is  within  890  km  of  the  true  position  and  the  receiver  is 
near  the  earth,  then  the  error  made  in  estimating  the  direction  cosines 
of  the  line-of-sight  unit  vector  Nj  will  then  be  less  than  10""^.  In  a 
second  approximation,  estimates  of  the  Doppler  correction  factor  which 
use  velocities  and  a  value  of  Nj  accurate  within  a  few  percent  should  be 
adequate. 

Estimating  receiver  velocity.  In  the  present  example  we  consider 
a  receiver  which  has  an  estimated  ground  velocity  of  Vg  =  300  m/s  in 
a  direction  30°  north  of  east,  parallel  to  earth's  surface.  Suppose  the 
receiver  at  the  instant  tc  —  240  s  exactly  is  at  the  position 


X, 


i?  COS  35° 

5  224  663.389  m 

iyG584  = 

0 
i?  sin  35° 

— 

0 
3  658  348.690  m 

which  has  been  used  in  our  earlier  examples.  The  net  velocity  will  be 
a  linear  combination  of  the  ground  velocity  and  velocity  of  the  point 
XivG584  due  to  earth  rotation.  From  geometric  considerations,  the  net 
receiver  velocity  at  time  tc  in  this  example  will  be 


vc 


-Vy  sin  30°  sin  35° 

■t;^  cos  30°  +  nei?cos35° 

Vg  sin  30°  cos  35° 


If  the  position  and  ground  velocity  were  accurately  known,  the  net  re- 
ceiver velocity  in  the  ECIF  would  then  be 

-86.036  465  453" 

Vc  =       640.796  091481      m/s. 

122.872  806  643 


(25) 


(26) 


(2- 
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Therefore  in  the  present  example,  we  take  the  velocity  estimate  in  the 
ECIF  to  be 

"-86" 
vc  =       641      m/s. 
123  J 

Estimating  receiver  position.  We  illustrate  the  calculation  of  the 
first-order  Doppler  correction  with  satellite  1  in  Table  7.  The  unit  vector 
N  is  computed  from  the  satellite  and  receiver  positions,  as  in  eq  (D.3): 


N  = 


\rR{tR)  -ftI^t)!  ' 


where  rR{tji)  and  ftI^t)  are  the  ECIF  position  vectors  of  receiver  and 
transmitter  at  reception  and  transmission  times  tR  and  tx-  At  the  re- 
ceive time,  with  bias  of  correction  423  ns,  the  receiver  position  vector  is 
actually 


(28) 


{D.3) 


rRi  = 


rc+vc(240.069  124  140  3  -  0.000  000  423  0  -  tc) 
5  224  656.660 


50.117 
3  658  358.300 


m 


(29) 


and  we  therefore  take  as  an  estimate  of  this  position  the  vector 


ri?i 


6  000  000 

23  000 

3  100  000 


m, 


(30) 


and  the  receiver  velocity  is  given  in  eq  (28),  which  was  discussed  above. 
The  bias  correction  was  not  significant  in  arriving  at  this  rough  estimate 
of  receiver  position. 

Estimating  satellite  velocity.  The  satellite  position  at  transmit  time 
tc  is  given  in  the  first  row  of  Table  7.  The  satellite  velocity  is  also  needed. 
It  may  be  obtained  by  differentiating  the  satellite  positions,  which  are 
derivable  from  the  prescription  of  Table  20-IV  in  the  ICD-GPS-200,  with 
respect  to  time,  and  then  applying  a  rotation  into  the  ECIF.  For  the 
present  example  the  satellite  velocity  in  the  ECIF  is 


Vi  = 


-304.267  760 
3  853.051  168 
-257.220  637 


m/s 


(31) 
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The  position  estimates  result  in  the  following  calculated  values  of  the 
components  of  Ni: 


Ni 


-0.437  588 
-0.114  125 
-0.891  904 


while  the  "true"  value  of  Ni  obtained  from  the  simulation  is 


(32) 


Ni 


-0.477  846  083 
-0.115  974  357 
-0.870  754  311 


(33) 


For  the  present  case,  Atsv  =  -0.000  056  554  783  4  s.  Then  eq  (16) 
gives  the  estimate 

Atm  ~  -0.000  056  554  784   s . 

This  estimate  agrees  extremely  well  with  the  "true"  value  of  Atjn  ob- 
tained from  the  simulation.  For  details  of  this  comparison,  see  Appendix 
E. 

Thus  the  Doppler  correction  factor  given  in  eq  (16)  can  be  computed 
by  the  receiver  and  added  to  tjij  to  predict,  with  excellent  accuracy,  the 
arrival  time  tjij  +  dtjij  of  the  code  which  corresponds  to  transmission 
of  tc  from  the  satellite.  Thus  only  times  of  arrival  of  code  boundaries 
would  actually  be  measured  by  the  receiver. 

Three  additional  calculations  using  eq  (16)  are  given  in  the  third 
column  of  Table  7  for  the  satellites  used  in  the  present  example.  Here 
the  time  tc  equals  240.000  000  000  s  exactly.  The  first  column  of  Table  7 
gives  local  clock  arrival  times  of  the  phase  time  r  =  tc-  The  second  col- 
umn lists  the  values  of  the  SV  clock  and  relativistic  corrections  derived 
from  the  navigation  messages.  The  third  column  lists  the  corrections 
AtRj  obtained  from  eq  (16).  Note  that  the  largest  difference  between 
Atsvj  and  Atnj  is  little  more  than  0.5  ns.  The  maximum  possible  such 
difference  is  about  4  ns  and  could  occur  for  a  clock  in  a  satellite  in  low 
earth  orbit.  The  resulting  values  of  local  clock  arrival  times  tcj  are 
listed  in  the  last  column  of  Table  7  and  are  repeated  in  the  last  column 
of  Table  8.  There  may  still  be  a  local  clock  bias  in  these  data.  In  this 
example,  the  receiver  clock  has  a  bias  of  —423  ns. 
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Table  7.  Computation  of  first-order  Doppler  corrections. 

(units  =  s) 


ID  #          tRj  AtsVj  ^tRj   eq  (16)  tcj   =  tRj   +  AtRj 

1  240.069  180  695  048  -0.000  056  554  796  2  -0.000  056  554  784  240.069  124  140  263 

2  240.078  615  674  081  0.000  493  792  700  6  0.000  493  792  194  240.079  109  466  275 

3  240.075  360  518  523  0.000  259  725  929  9  0.000  259  725  401  240.075  620  243  923 

4  240.072  939  857  576  0.000  305  714  376  7  0.000  305  713  988  240.073  245  571  564 


Table  8.  Satellite  positions  r^  and  mecisiured  arrival  times  tcj- 


I 

15  126  951.488 

2 

20  996  573.361 

3 

23  391  384.147 

4 

10  656  983.742 

ID  #           Satellite  Xj  (in)            Satellite  IJj  (m)  Satellite  Zj  (m)  Local  clock  time  tcj  (s) 

2  403  354.115  21  702  797.709               240.069  124  140  3 

-15  349  939.783  -5  178  329.880               240.079  109  466  3 

10  956  612.394  —4  333  140.777                240.075  620  243  9 

-11  270  652.822  21  703  578.299               240.073  245  571  6 


Table  8  gives  the  data  from  the  all  four  measurements:  the  derived 
positions  of  the  satellites  at  GPS  system  transmission  times  and  the 
corresponding  arrival  times  tcj  on  the  receiver  clock. 

(6)  The  iterative  solution,  which  has  been  carried  out  as  described  in  Ap- 
pendix C  is  given  at  successive  stages  in  Table  9.  At  the  initial  stage,  estimates 
for  the  ECEF  receiver  position  are  inserted.  For  the  sake  of  illustration,  we  as- 
sume the  initial  estimate  of  the  receiver  position  places  it  at  earth's  center.  The 
value  of  the  local  clock  bias  is  initially  estimated  to  be  300  ns.  The  iteration  pro- 
cess converges  rapidly,  to  within  a  millimeter  or  so  of  the  correct  ECEF  receiver 
position. 

(7)  No  rotations  are  necessary  in  the  final  stage  since  the  position  has  been 
determined  in  the  ECEF  frame. 

The  receiver  positions  at  the  times  tcj  of  reception  might  be  needed  for  other  purposes. 
These  could  be  obtained  using  a  best  estimate  of  the  velocity  of  the  receiver  in  the  ECIF 
(the  ECEF  frame  at  tc)  using  eq  (18).  The  components  of  re  are  in  this  reference  frame. 
A  second  series  of  measurements  at  a  slightly  later  time  t^  will  give  a  new  value  of  tq  in 
the  new  reference  frame.  These  two  values  could  be  differenced  or  filtered  to  provide  an 
update  to  the  ground  velocity.  Such  an  updated  value  could  be  transformed  into  an  ECIF 
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Table  9.  Iteration  solution  for  re  and  receiver  clock  bias  be- 


Trial  # 

Receiver  Xq  (m) 

Receiver  yc   (m) 

Receiver  Zc   (m) 

Local  clock  bicus  be  (s) 

0 

0.0 

0.0 

0.0 

0.000  000  300  0 

1 

5  051  460.568 

-9  835.265 

3  545  131.610 

-0.002  453  425  2 

2 

5  226  991.162 

132.182 

3  659  870.280 

0.000  032  544  3 

3 

5  224  663.799 

-0.023 

3  658  349.573 

-0.000  000  417  2 

4 

5  224  663.390 

-0.000 

3  658  348.690 

-0.000  000  423  0 

5 

5  224  663.390 

-0.000 

3  658  348.690 

-0.000  000  423  0 

with  the  help  of  a  rotation  such  as  in  eq  (28),  to  be  used  in  further  position  measurements. 

Example  4.52    Time  Tagging  at  Transmitter  Neglecting  Doppler  Corrections 

The  first-order  Doppler  corrections  to  Atsv  are  very  small,  as  seen  in  Example  4.51, 
above.  Here  we  repeat  Example  4.51  with  one  change:  the  Doppler  corrections  are  com- 
pletely neglected.  Since  including  these  first-order  Doppler  corrections  usually  amount  to 
changing  the  arrival  times  on  the  local  clock  by  less  than  1  ns,  only  a  few  centimeters  of 
error  would  be  expected  to  occur  when  such  corrections  are  neglected.  Table  10  gives  the 
satellite  positions  Vj  and  the  measured  arrival  times  tcj  when  this  approximation  is  made. 
The  table  is  identical  to  Table  8  except  in  the  fourth  column. 

Solutions  of  Eqs.  (21)  by  iteration  are  given  in  Table  11.  Thus  neglecting  the  first- 
order  Doppler  correction  causes  a  position  error  of  about  75  cm  in  this  case. 

Table  10.  Satellite  positions  Vj  and  arrival  times  tcj  neglecting  Doppler  correction. 

ID  #  satellite  Xj  (m)  satellite  yj  (m)  satellite  Zj  (m)  tcj  ~  tjij  -\-  /\tsvj  (s) 


1 

15  126  951.488 

2  403  354.115 

21  702  797.709 

240.069  124  140  3 

2 

20  996  573.361 

-15  349  939.783 

-5  178  329.880 

240.079  109  466  8 

3 

23  391  384.147 

10  956  612.394 

-4  333  140.777 

240.075  620  244  5 

4 

10  656  983.742 

—  11  270  652.822 

21  703  578.299 

240.073  245  572  0 

5.  TIME  TRANSFER  RECEIVER 

In  this  case  we  may  use  each  pseudorange  measurement  from  each  satellite  sepai^ately 
to  estimate  our  clock  offset  from  GPS  time.  We  either  time  tag  measurements  at  a  trans- 
mission time  and  measure  time  of  reception,  or  time  tag  measurements  at  a  reception  time 
and  determine  time  of  transmission  from  the  code  lock  of  the  receiver.   We  then  use  the 
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Table  11.  Solution  for  re  and  receiver  clock  bias  be  neglecting  Doppler  correction. 


Trial  # 

Receiver  Xc  (m) 

Receiver  yc  (m) 

Receiver  Zq  (m) 

Local  clock  bias  be  (s) 

0 

0.0 

0.0 

0.0 

0.000  000  300  0 

1 

5  051  461.207 

-9  835.359 

3  545  131.862 

-0.002  453  427  6 

3 

5  226  991.866 

132.086 

3  659  870.566 

0.000  032  542  4 

4 

5  224  664.501 

-0.072 

3  658  349.243 

-0.000  000  419  1 

5 

5  224  664.092 

-0.096 

3  658  348.976 

-0.000  000  424  9 

6 

5  224  664.092 

-0.096 

3  658  348.976 

-0.000  000  424.9 

prescriptions  from  the  ICD-200,  Table  20-IV,  to  obtain  the  satellite  position  in  the  ECEF 
frame.  We  must  estimate  the  geometric  path  delay  t^  in  addition  to  the  ionospheric 
and  tropospheric  corrections.  Since  we  know  the  receiver  position  we  may  compute  t£)  as 
follows. 

If  ECI  position  vectors  are  referenced  to  the  time  of  signal  transmission,  then 

\v{tn)-Ii{tT)\  ^Jr{tT)-mT)\    ,    [r{tT)-n{tT)]-v  ._ 

c  c  c^ 

where 

R  is  the  SV  position  vector, 

r  is  the  receiver  position  vector, 
tT  and  tn  are  GPS  time  at  respectively  transmit  and  receive  times,  and 

V  is  the  receiver  velocity  vector  in  the  ECIF. 

The  last  term  in  eq  (34)  is  sometimes  called  the  Sagnac  correction.  In  the  case  of  an 
earth-fixed  user, 

V  =  i7e  X  r.  (35) 

If  ECIF  position  vectors  are  referenced  to  the  time  of  signal  reception  (instead  of 
satellite  transmission  time),  then 

\v{tR)-K{tT)\  ^^  \r{tR)-K{tn)\       [r(t;,)  -  R(^h)]  •  V         .... 

to  =tR  —  tT  =  « 5 (3b) 

c  c  c^ 

where: 

V  is  the  SV  velocity  vector  in  the  ECI. 

A  single- channel  receiver  with  antenna  placed  at  a  known  location  can  provide  a  means 
for  accurate  transfer  of  GPS  time  to  a  local  clock.  Let  the  antenna  be  at  the  position 
^ECEF  whose  earth-fixed  coordinates  are  accurately  known  and  do  not  change  with  time. 
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The  appropriate  expression  for  the  total  propagation  delay  for  a  signal  from  one  satellite 
is  eq  (34).  In  eq  (34)  the  positions  are  ECI  positions  but  are  referred  to  the  time  of 
transmission  t^.  The  first  term  on  the  right  side  of  eq  (34)  is  proportional  to  the  distance 
between  receiver  and  transmitter  at  the  transmit  time  tx',  if  we  call  this  distance  L,  then 

L  =  \r{tT)-Ti{tT)\.  (37) 

Although  the  position  vectors  r{t'r)  and  R(tT)  are  ECI  position  vectors,  the  distance 
L  is  invariant  under  spatial  rotation  and  can  therefore  be  computed  in  any  convenient 
coordinate  system.  The  ECEF  position  of  the  satellite  at  transmit  time  is  directly  derivable 
from  the  prescription  given  in  the  ICD-200,  Table  20-IV,  and  therefore  no  rotations  are 
necessary  to  compute  L: 

L  =  IrECEFitr)  —'R'ECEF{tT)\-  (38) 

The  second  term  in  eq  (34),  the  Sagnac  correction  term,  can  be  further  reduced  for  this 
case  since  it  also  is  a  scalar  with  respect  to  rotations.  The  ECI  velocity  of  the  earth-fixed 
receiver  is  due  to  earth  rotation  and  is 

V  =  17e  X  r{tT).  (39) 

The  Sagnac  correction  term,  which  we  denote  by  Atsagnaci  is  therefore 

[r{tT)-R{tT)]-^e^r{tT) 


At 


Sagnac 


tie 

■r{tT)  X 

C2 

[r{tT)-K{tT)] 

^e 

■r{tT)x 

C2 
[-R(iT)] 

He 

C2 

•  K{tT)  X 

T{tT) 

(40) 
X  r{tT) 

c2 

^e  •  ^ECEPitT)  X  ^ECEFitr) 


c2 


l^el  {X EC E F{iT)y EC EFitT)  —  yECEF{tT)x ECEFitr)) 


C2 


where  standard  vector  operations  have  been  used  to  simplify  the  scalai*  product.  In  the 
last  step,  the  only  component  of  the  cross  product  which  contributes  is  parallel  to  the 
earth's  angular  velocity  vector;  this  product  is  invariant  under  rotations  about  this  axis, 
so  the  computation  can  be  performed  in  the  ECEF  frame  at  the  instant  of  transmission. 
Thus  for  a  timing  receiver  no  rotations  need  be  applied  to  compute  the  local  time. 
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The  final  result  in  eq  (40)  above  can  be  put  in  a  form  which  is  easy  to  interpret  and 
can  be  found  in  many  places  in  the  literature: 

AtSagnac  =  "^^S  ,  (41) 

where  Ae  —  |  [R(iT)  x  ^{tT)]z  i^  the  equatorial  projection  of  the  area  swept  out  by  a 
vector  from  the  center  of  the  earth  to  the  path  of  the  timing  pulse  as  it  propagates  from 
transmitter  at  R(iT)  to  receiver  at  r{tT)-  Clearly  this  area  can  be  calculated  in  either 
ECIF  or  ECEF  coordinates. 

Example  5.1    Earth-fixed  Single- Channel  Timing  Receiver 

We  give  four  examples  here  based  on  the  data  from  four  different  satellites  given  in 
Table  1.  For  all  these  examples  the  ECEF  position  of  the  receiver  is  that  specified  in  eq 

(7). 

In  Table  12  the  first  column  labels  different  satellites  corresponding  to  the  data  entries 
in  Table  1.  The  second  column  gives  the  transmission  epochs.  The  third  column  gives  the 
transmission  delay  contributions  X/c,  the  fourth  column  gives  the  Sagnac  corrections, 
and  the  fifth  column  sums  the  entries  in  the  previous  three  columns  to  give  the  time  at 
the  receiver.  The  agreement  of  the  entries  in  the  last  column  refiects  consistency  of  the 
procedure  across  the  satellite  constellation. 

Table  12.  Examples  of  time  transfer  to  a  single-channel  earth-fixed  receiver. 

SV  #    Transmission  epoch  ij-  (s)      L/C  (s)        ^^Sagnac  (s)  tjij   (s) 

1  37  239.924  422  365  6  0.075  577  714  8  -0.000  000  080  5  37  240.000  000  000  0 

2  37  239.920  713  391  8  0.079  286  677  5  -0.000  000  085  5  37  240.000  000  000  0 

3  37  239.925  307  870  0  0.074  692  197  4  -0.000  000  067  4  37  240.000  000  000  0 

4  37  239.929  346  353  9  0.070  653  609  2  0.000  000  036  9  37  240.000  000  000  0 
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Appendix  A. 
The  Eccentricity  Correction 

The  epoch  of  the  transmission  event  from  each  satelhte  must  be  corrected  to  account 
for  variations  in  the  proper  time  of  each  satellite  clock  due  to  orbital  eccentricity.  This 
aspect  of  system  design  of  the  GPS  was  fixed  very  early  when  satellites  did  not  have  much 
computing  power.  The  correction  compensates  the  atomic  clock  in  the  satellite  for  periodic 
gravitational  frequency  shifts  and  second-order  Doppler  shifts  due  to  the  satellites'  periodic 
change  of  radius  and  velocity  in  its  eccentric  orbit.  The  correction  converts  the  proper  time 
kept  on  the  atomic  clock,  to  coordinate  time — that  is,  to  GPS  time.  To  put  it  another 
way,  a  satellite  in  an  eccentric  orbit  samples  different  values  of  earth's  gravitational  field 
as  the  satellite  moves  toward  and  away  from  the  earth.  Consequently  it  suffers  a  periodic 
gravitational  frequency  shift.  At  the  same  time,  conservation  of  angular  momentum  causes 
the  satellite  to  speed  up  or  slow  down,  giving  rise  to  a  periodic  second-order  Doppler  shift. 
Thus  SV  clocks  must  be  corrected  for  these  effects  so  that  the  corrected  clocks  effectively 
run  at  a  constant  rate,  the  rate  of  coordinate  time. 

It  is  important  to  realize  that  this  correction  is  the  same  for  all  receivers,  no  matter 
where  they  are  or  what  their  state  of  motion  is.  The  correction  relates  solely  to  the 
transmitter.  In  particular  it  is  correct  for  users  of  GPS  anywhere  in  space.  It  is  incorrect 
to  claim,  as  is  done  in  Ref.  [6],  that  the  correction  applies  only  to  users  near  earth's 
surface. 

Let  us  denote  the  time  that  the  user  estimates  from  the  information  encoded  in  the 
signal  as  tsv  a^  in  the  ICD-GPS-200.  This  must  be  corrected  according  to  the  section 
describing  the  user  algorithm  to  obtain  GPS  system  time  t  (Sect.  20.3.3.3.3.1  in  the  ICD- 
GPS-200C).  The  correction  has  several  different,  but  equivalent  forms.  Using  the  notation 
of  ICD-GPS-200C,  Sect.  20.3.3.3.3.1,  the  relativistic  correction  to  the  transmission  epoch 

tsv  is 

^VGMA     .     ^ 

Atr  =  -  2 ^ esmEk 

c^ 

=  Fe^sinEk  (^-1) 

2R   V 

where 

F  =  -4.442  807  633  X  10"^°  s//l^,  (^-2) 

and  where  GM  is  the  product  of  the  Newtonian  gravitational  constant  and  the  mass  of 
the  earth,  A  is  the  semi-major  axis  of  the  satellite  in  meters,  e  is  the  eccentricity,  Ek 
is  the  satellite's  eccentric  anomaly,  R  is  the  instantaneous  position  vector  of  the  SV  at 
transmission  time,  and  V  is  the  instantaneous  velocity  vector  of  the  SV.  We  need  GPS 
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time  t  before  we  can  obtain  Ek  as  prescribed  in  Table  20- IV,  so  the  equations  are  coupled. 
It  suffices  to  use  tsv  to  determine  E^,  then  obtain  GPS  time  t.  Also,  since  the  scalar 
product  in  the  last  form  of  the  correction  is  invariant  under  spatial  rotations,  it  can  be 
evaluated  in  any  inertial  frame.  This  correction  Atr  is  subtractedhom  the  transmitted  time 
tsv  to  obtain  GPS  system  time  t.  The  time  t  may  then  be  applied  to  the  user  algorithm 
for  ephemeris  (Table  20-IV  in  the  ICD-GPS-200C),  to  obtain  the  x,  y,  z  coordinates  of  each 
satellite  at  transmission  time. 

Derivations  of  the  above  results  can  be  found  in  references  [2],  [3],  and  [4]. 
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Appendix  B. 

Linearization  Algorithm  for 
Time  Tagging  at  Receiver 

We  now  explain  one  possible  method  of  linearization  of  the  one-way  light  time  equa- 
tions which  are  valid  in  an  inertial  coordinate  system. 

Consider  four  GPS  satellites  which  at  GPS  times  ti,  ^2,  ^3,  ^4  send  out  signals  from 
the  ECIF  locations  ri,  r2,  r3,  r^.  We  assume  that  these  data  have  been  expressed  in  the 
chosen  ECI  reference  frame.  These  four  signals  are  received  simultaneously  at  time  t  by  a 
receiver  at  position  r.  The  problem  is  to  determine  t  and  r  at  the  receiver.  The  velocity 
of  the  receiver  does  not  have  to  be  considered.  The  receiver's  position  r  at  the  time  of 
reception  t  will  be  determined  by  the  solution  of  four  equations  which  express  the  condition 
that  the  speed  of  propagation  is  c.  The  four  equations  are 

{r-rjf-c\t-t-,f=0;     i  =  1,2,3,4.  (B.l) 

Such  equations  are  clumsy  to  solve  in  general.   However,  if  the  position  r  and  time  t  are 
known  approximately,  the  equations  can  be  reduced  to  a  system  of  four  linear  equations 
which  can  be  solved  by  standard  matrix  inversion  techniques. 
Suppose  that  for  user  position  and  time  we  write 

r  =  r(*)+Ar;      t  =  t^'^  +  At,  {B.2) 

where  the  quantities  r^*^,  t^'^^  represent  guesses  or  trial  values  for  the  position  and  time, 
and  Ar,  At  are  assumed  to  be  small  quantities  which  are  to  be  determined  by  solving  eqs 
(B.l).  The  index  i  in  eq  (B.2)  labels  the  zth  trial  value  in  a  repeated  or  iterative  process 
of  solving  Eqs.  (B.l). 

We  substitute  eqs  (B.2)  into  eqs  (B.l),  expand  and  neglect  squares  of  the  small 
quantities  Ar  and  At.  Then  the  equations  take  the  following  form 

(r(*)  -  r,)  •  Ar  -  ^{t^'^  -  tj)At  =  -  ((P{t^'^  -  tjf  -  (r^''  -  r,)^)  ,  {B.3) 

for  j  —  1,  2,  3, 4.  In  eqs  (B.3),  the  quantity  c(t^*^  ~^j)  is  the  ith  estimate  of  the  pseudorange 
from  the  receiver  to  the  jfth  satellite.  Eqs  (B.3)  are  a  system  of  linear  inhomogeneous 
equations  which  can  be  solved  by  matrix  inversion.  Numerical  examples  of  the  solution 
are  discussed  in  the  main  text.  The  matrix  of  coefficients  of  the  unknowns  Ar  and  At 
will  usually  be  nonsingular,  unless  the  configuration  of  satellites  is  so  unfavorable  that  eqs 
(B.l)  do  not  have  a  solution.  The  solutions  obtained  will  be  approximate,  but  can  be  used 
to  obtain  new  trial  values;  the  process  of  iteration  can  then  be  repeated. 
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Thus,  suppose  at  stage  i  in  the  iteration,  eqs  (B.3)  have  been  solved  for  corrections 
Ar  and  At;  then  better  approximations  can  be  found  by  setting 

r(*+i)  ^  j.{r)  _^  ^j.;      i(^+i)  =  t(^)  +  At.  (BA) 

The  iteration  converges  rapidly,  particularly  if  the  initial  trial  guesses  were  close.  The 
solutions  will  be  found  at  some  stage  when  the  inhomogeneous  terms  on  the  right  sides  of 
eqs  (B.3)  become  negligible.  Then  the  solutions  satisfy 

(rW_r^.)2_c2(i(0_t^,)2^0;     i  =  1,2,3,4,  (R5) 

which  are  equivalent  to  eqs  (B.l). 

In  summary  the  computation  scheme  proceeds  as  follows: 

Step  1.  Choose  approximate  values  for  user  position  and  time. 

Step  2.  Solve  eqs  (B.3)  to  obtain  corrections  Ar  and  At  for  user  position 
and  time. 

Step  3.  Test  whether  the  corrections  are  sufficiently  small  that  computa- 
tion can  be  halted;  if  not,  combine  the  corrections  obtained  in  Step  2  with  the 
approximate  values  of  position  to  obtain  a  new  approximate  position.  Then  go 
back  to  Step  2. 
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Appendix  C. 

Linearization  Algorithm  for 
Time  Tagging  at  Transmitters 

When  the  measurements  are  time  tagged  at  the  transmitters,  then  the  measurements 
at  the  receiver  will  occur  at  unequal  times.  This  entails  a  different  method  for  linearization 
of  the  propagation  delay  equations  than  that  discussed  in  Appendix  B. 

Suppose  that  we  decide  to  time  tag  the  signals  from  four  GPS  satellite  transmitters 
at  the  common  GPS  system  time  tc-  Let  the  arrival  times  of  these  signals  on  the  receiver 
clock  be  tjij,  and  let  the  positions  of  the  sateUites  at  time  tc  from  which  these  signals 
were  transmitted  be  r^,  with  j  =  1, 2, 3, 4.  These  positions  are  derived  from  the  navigation 
messages;  they  will  be  expressed  in  the  ECEF  frame  at  the  instant  tc-  Therefore  it  is 
highly  advantageous  to  work  with  an  ECI  reference  frame  which  coincides  with  the  ECEF 
frame  as  though  it  were  frozen  at  tc-  Then  no  rotations  will  be  necessary  to  bring  the 
satellite  positions  into  a  common  ECI  frame — they  are  already  expressed  in  the  ECIF 
defined  at  tc-  If  we  solve  the  propagation  delay  equations  for  the  receiver  position  at  tc 
we  have  the  receiver  position  in  this  same  ECIF.  Thus  time  tagging  at  the  transmitters 
.has  the  great  advantage  that  no  rotations  are  needed. 

In  general  the  signal  propagation  times  from  satellites  to  receiver  will  not  be  equal — 
they  can  vary  roughly  from  67  to  86  ms  for  a  receiver  on  earth's  surface.  This  means  that 
it  will  be  necessary  to  have  an  estimate  of  the  receiver  velocity  at  the  time  tc  in  order  that 
the  different  receiver  positions  at  the  different  signal  arrival  times  can  be  accounted  for. 
We  may  also  need  the  receiver's  acceleration.  We  assume  these  are  available,  to  sufficient 
accuracy,  from  previous  position  solutions  or  from  a  filtering  algorithm. 

We  will  illustrate  this  here  by  assuming  the  receiver  position  is  expressible  as  a  Taylor 
series 

r/?(t)-rc  +  vc(i-ic)  +  -,  (C.l) 

where  re  and  vc  are  the  receiver  position  and  velocity  in  the  chosen  ECIF,  at  time  tc-  An 
additional  acceleration  term  in  this  expansion  might  be  needed  for  a  receiver  in  a  highly 
dynamic  trajectory.  We  assume  here  that  vc  is  known  but  r^  is  unknown. 

The  measured  arrival  times  tjij  of  the  signals  on  the  receiver  clock  will  in  general  not 
be  true  GPS  system  times  because  of  a  possible  bias  be  on  the  receiver  clock.  Thus  the 
GPS  system  times  of  arrival  will  be 

tj  =  tRj  +  be.  (C.2) 

and  as  part  of  the  data  analysis,  the  bias  will  be  found.  The  receiver  positions  at  times  tj 
will  therefore  be 

mitj)  =  TRj  =rc  +  ycit'Rj  +bc-  tc)  ■  {C.2) 
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In  this  Appendix  we  do  not  discuss  code  boundary  measurements  and  corrections.   This 
important  topic  is  discussed  in  detail  in  Sect.  3.  We  simply  assume  that  the  appropriate 
places  in  the  code  have  been  identified  by  the  receiver. 
The  four  propagation  delay  equations  are  then 


(VRj  -Tjf  ^  C^itRj  +bc  -tcf 


(C.4) 


or 


(re  +  y^citRj  +bc-  tc)  -  Yjf  =  c^{tRj  +  bc-  tcf  . 


(C.5) 


The  unknown  position  vc  and  bias  be  may  be  found  from  these  four  equations  cis  follows. 
As  in  Appendix  B,  we  assume  we  can  write  the  unknown  quantities  in  terms  of  good 
estimates,  plus  small  corrections.  Thus  we  write 


re 


c 


+  Ar,      be  =  b^fi  +  Ai , 


(C.6) 


where  the  quantities  Yq\  b^  represent  trial  estimates  for  the  position  vector  and  clock 
bias,  and  Ar,  At  are  small  quantities.  The  index  (i)  in  eq  (C6)  labels  the  i*'*  trial  value 
in  a  repeated  or  iterative  process  of  solving  eqs  (C.5). 

We  substitute  {C.Q)  into  (C.5),  expand  and  neglect  squares  of  small  quantities  Ar 
and  At.  Then  the  equations  reduce  to 


R 


{^) 


Ar-c' 


^(.)  _  ^c  ■  R; 


At  = 


1 


T. 


(') 


-    R 


M) 


(C.7) 


where 


and 


Tf^  =  tRj  +  6^'  -  tc  , 


Rf  =  r^;)  +  ycT^;^ 


(C.8) 
(C.9) 


This  system  of  equations  is  of  the  same  form  as  that  discussed  in  Appendix  B  and  can  be 
solved  by  iteration  in  a  similar  manner  as  that  discussed  at  the  end  of  Appendix  B.  The 
main  difference  from  eq  (B.3)  is  that  the  velocity  estimate  vc  appears  in  eq  (C.7)  in  three 

(i) 

places:  both  on  the  left  and  the  right  sides  in  the  definition  of  R^  ,  and  in  the  coefficient 
of  At  on  the  left  side. 

Thus  the  iteration  scheme  proceeds  as  follows. 

Step  1.  Choose  tc;  estimate  vq  in  the  ECIF.  This  will  be  a  combination  of 
velocity  relative  to  the  ground,  and  velocity  of  the  ground  relative  to  the  chosen 
ECIF,  due  to  earth  rotation. 
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Step  2.  Choose  approximate  values  for  the  user  position  re  and  clock  bias 
he.  Compute  ¥if  and  T^/K 

Step  3.  Solve  eqs  (C7)  to  obtain  corrections  Ar  and  At  for  user  position 
and  clock  bias. 

Step  4.  Test  whether  the  corrections  are  sufficiently  small  that  computa- 
tion can  be  halted;  if  not  combine  the  corrections  obtained  in  Step  3  with  the 
approximate  values  of  position  and  bias,  to  obtain  a  new  approximate  position 
and  bias. 

Step  5.  Recompute  R^*^  and  T-*^  and  return  to  Step  3. 
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Appendix  D. 

Equivalence  of  Integrated  Doppler 
and  One- Way  Light  Time 

In  this  Appendix  we  show  the  complete  equivalence  of  the  one-way  light  time  mea- 
surement to  integrated  Doppler  measurements. 

Neglecting  gravitational  path  bending  and  gravitational  time  delay  is  well  justified; 
electromagnetic  signals  travel  with  speed  c  in  straight  lines  in  the  ECI  reference  frame. 
The  direction  of  propagation  may  then  be  described  by  a  constant  unit  vector  N  along  the 
path  from  transmitter  to  receiver. 

Let  the  coordinate  time  of  transmission  of  a  signal  be  denoted  by  tT  and  the  coordi- 
nate time  of  reception  be  denoted  by  tji.  The  positions  of  the  transmission  and  reception 
events  in  ECI  coordinates  are  denoted  by  rT{tT)  and  rR{tR).  Also,  the  velocities  of  the 
transmitter  and  receiver  at  the  transmission  and  reception  epochs  are  denoted  by  V  and  v, 
respectively;  this  notation  is  introduced  in  eqs  (1-2)  where  the  Doppler  effect  was  stated. 
In  this  Appendix  we  derive  the  Doppler  formula  including  gravitational  frequency  shifts. 

The  one-way  propagation  of  a  signal  from  transmitter  to  receiver  gives 

tR  =  K  +  tT  +  -IrnitR)  -  ftI^t)!  ,  {D.l) 

where  K  represents  a  possible  systematic  error  or  bias  in  the  coordinate  time  of  the  receiver 
clock.  The  term  \rji{tfi)  —  rT{tT)\/c  is  the  propagation  delay  in  the  equation,  which  must 
be  accounted  for  in  any  use  of  GPS  signals  for  navigation  or  time  transfer.  We  show  here 
that  if  a  receiver  properly  accounts  for  this  term,  then  it  has  automatically  accounted  for 
the  relativistic  Doppler  effect,  integrated  over  time. 

Let  there  be  a  second  signal  sent  from  transmitter  to  receiver  after  one  period  of 
oscillation  of  the  oscillator  in  the  transmitter.  Let  the  times  of  transmission  and  reception 
of  this  second  signal  be  tT  +  dtx  and  tn  -\-  dtn.  The  frequencies  will  be  very  high,  so  the 
small  differentials  may  be  calculated  by  taking  differentials  of  eq  (D.l).  The  term  in  K 
drops  out,  and 

dtR  =  dtT+  --. — 7-— -—-  -{drR-drT)  ,  {D.2) 

c\rR{tR)  -rritTn 

where  dvR  and  drx  are  the  changes  in  positions  of  the  receiver  and  transmitter  during  the 

time  intervals  dtR  and  dix-  But  the  unit  vector  N  is  exactly 

\rR{tR)  -rT{tT)\ 
and  during  the  small  time  intervals  dtR,  dtx  we  have 

dvR  =  \dtR,      dvT  =  VdtT  .  {DA) 

42 


Therefore 

dtR  =  dtr  +  -N  •  {\dtR  -  VdtT)  •  •  (D.5) 

Rearranging  gives 

,      /        N-v\        ,      /        N-V\ 
dtR  (  1 —  j  =  dtT  (  1 —j  .  {D.6) 

Next  we  relate  the  coordinate  time  intervals  dtR,  dtT  to  the  proper  frequencies  of 
received  signal  and  transmitted  signal.  This  part  of  the  calculation  is  only  carried  to  terms 
of  order  c~^ .  Let  $(r)  be  the  Newtonian  gravitational  potential  at  position  r  in  ECI 
coordinates.  Let 

A$  =  $-$o,  (^-7) 

where  the  constant  $o  is  the  effective  gravitational  potential  on  the  rotating  geoid,  in- 
cluding the  centripetal  contribution  arising  from  earth's  spin.  Then  the  proper  time  dr 
that  elapses  on  an  atomic  clock  as  it  moves  through  the  position  increment  dr  during  the 
coordinate  time  increment  dt  is  given  approximately  according  to  general  relativity  by  [2,3] 


c^dT^  =[l  +  t^\c^dt^-  {dvf  =  h  +  ^  -  i^^^  )  c^dt^  .  {D.8) 


This  may  be  solved  for  the  coordinate  time  increment  dt  in  terms  of  the  proper  time: 

dr  dr  .  _     . 

dt  =      ,  ~  — rr-JiTTT  ■  {D.9) 

v/(l  +  2A$/c2-(drM)2/c2)        l  +  ^-liJLMll 

where  the  square  root  has  been  expanded  using  the  first  term  in  a  binomial  expansion, 
y/l  +  x/c^  ~  1  +  x/{2c^).  This  may  be  applied  to  the  coordinate  time  intervals  dtR  and 
dtx  between  receptions  and  transmissions  of  the  two  successive  signals.  Accurate  to  terms 
of  order  1/c^, 

dtR^  rx7-T T-  (D.IO) 

dTT 
dtT  =  T^n-A T-  (^11) 

■^  "•"         c2  2    c2 

We  substitute  these  results  into  eq  (D.6)  and  obtain 

dra  (1  -  ^)  drr  (l  -  ^) 


.         A^(rfl)    _  IyI  1     I     A^(rT)   _   1  V^ 

■^    "•"  C2  2    C2  ^    '^  C^  2     C2 


{D.12) 


The  last  step  is  to  relate  the  proper  time  intervals  to  the  proper  frequencies.  Using 
F  and  /  for  the  proper  transmitted  and  received  frequencies,  respectively,  we  find 

{D.13) 


1 

dTT 

.        f  = 
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1 
dTR 

And  substitution  into  eq  (D.12)  gives  the  final  relation  between  proper  frequencies: 

/  (1  +  A$(ri?)/c2  -  v2/2c2)        F  (1  +  A$(rT)/c2  -  Y'^ I2c^) 


1_N:X  -  1_N^  •  (^-14) 


c 


Upon  neglecting  the  gravitational  frequency  shifts,  this  form  agrees  with  eq  (1)  in  the  main 
text  to  the  required  order  of  approximation. 

Several  observations  about  the  result  in  eq  (D.14)  are  appropriate.  First,  all  quanti- 
ties, except  the  proper  frequencies,  are  measured  in  ECI  coordinates.  Second,  there  are  no 
terms  of  the  form  v  •  V/c^.  Third,  to  obtain  the  ratio  f  I F  or  F//  correct  to  order  1/c^,  we 
must  take  care  to  expand  the  aberration  factors  (1  —  N  •  \/c^)~^  and  (1  —  N  •  V/c^)~^  to 
the  required  order.  To  illustrate  this  last  point,  let  us  calculate  the  ratio  F//.  Expanding 
eq  (D.14)  to  second  order  gives 

P/,       ,  ,^(rH)-^(rT)       v2        V2 

/  x2  (^-15) 

N-   v-V        N-vN-V       /N-v^^  ^         ^ 

+ ^ +  ' 


Conservation  of  coordinate  frequency.  Consider  the  special  case  when  both  the  trans- 
mitter and  the  receiver  are  at  rest  in  the  ECI  reference  frame.  Then  eq  (D.6)  gives 

dtji  =  dtT . 

This  means  that  the  coordinate  time  between  reception  of  successive  wavefronts  is  the 
same  as  the  coordinate  time  between  the  emission  of  successive  wavefronts.  The  coordinate 
frequency  is  the  reciprocal  of  this  time,  so  this  result  states  that  for  observers  at  rest  in  the 
ECI  reference  frame,  coordinate  frequency  is  conserved.  This  is  a  well-known  fact  that  can 
be  proved  by  separate  arguments  whenever  the  metric  tensor  coefficients  do  not  depend 
on  the  time.  This  is  the  case  here. 

We  may  readily  reverse  the  argument  given  above  and  pass  from  eq  (D.14)  back  to 
eq  (D.l)  by  integration  rather  than  differentiation.  Equation  (D.6)  is  an  intermediate  step 
in  such  an  argument,  which  is  obtained  by  combining  eqs  (D.IO,  D.ll,  and  D.12)  so  as 
to  eliminate  the  proper  time  intervals  drji  and  drT-  Then  eq  (D.2)  follows  from  (D.3-5). 
Integration  then  yields  eq  (D.l)  with  an  unknown  constant  of  integration. 

Thus  the  change  in  the  one-way  light  time,  as  used  in  GPS  for  pseudorange,  is  equal  to 
the  integrated  Doppler  frequency  shift  over  the  interval  of  the  change  using  the  frequency 
relation  (D.13). 
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Appendix  E. 

Numerical  Check  on 
the  First-Order  Doppler  Correction,  Equation  (16) 

In  this  Appendix  we  use  data  provided  in  Example  5.1  to  compare  the  "true"  value 
of  Atjij  from  the  simulation  with  values  obtained  using  eq  (16).  For  the  case  of  satellite  1 
in  Example  5.1, 

Atsv  =  -0.000  056  554  796  2  s. 

Then  eq  (16)  gives 

Atm  =  -0.000  056  554  784   s . 

This  estimate  using  first-order  Doppler  on  the  coordinate  time  may  be  compared  with 
the  simulation  which  gives  for  the  two  arrival  times  of  the  position  in  the  code  train 
corresponding  to  tsv  and  tc  the  following  time  difference: 

Atm  =240.069  124  140  263  -  240.069  180  695  048  s 
=  -  0.000  056  554  785  s . 

The  extremely  good  agreement  between  these  numbers — one  calculated  using  the 
Doppler  formula  and  one  using  propagation  delay  equations — suggests  the  following  strat- 
egy for  measuring  the  coordinate  time  of  arrival  of  the  code  corresponding  to  the  chosen 
measurement  time  tc-  If  we  choose  a  time  tc  as  some  multiple  of  1  ms  so  that  it  lies 
on  a  code  boundary,  the  SV  clock  correction  Atsvj  can  be  computed  from  the  naviga- 
tion message.  The  receiver  can  then  measure  the  arrival  time  tjij  of  the  code  boundary 
tsv  J  which  is  equal  numerically  to  tc  (arrival  of  phase  time  T  =  tc)-  Then  the  Doppler 
correction  factor  given  in  eq  (16)  can  be  computed  by  the  receiver  and  added  to  tRj  to 
predict,  with  excellent  accuracy,  the  arrival  time  tjij  +  dtjij  of  the  code  which  corresponds 
to  transmission  of  tc  from  the  satellite.  Thus  only  times  of  arrival  of  code  boundaries 
would  actually  be  measured  by  the  receiver. 

Three  additional  calculations  based  on  eq  (16)  are  given  in  Table  E.l,  for  the  four 
satellites  used  in  Example  5.1.  The  simulated  intervals  Atsvj  are  given  in  the  first  column 
of  Table  E.l.  The  second  column  lists  the  values  of  At^j  obtained  using  eq  (16).  For 
comparison,  the  last  column  gives  the  "true"  values  of  Atjij  obtained  from  the  simulation. 
In  all  cases  the  first-order  Doppler  correction,  eq  (16),  agrees  extremely  well  with  the 
"true"  values.  The  corrections  have  been  incorporated  into  the  data  listed  in  Table  7. 
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Table  E.l.  Computation  of  first-order  Doppler  corrections. 

(units  =  s) 


ID  #  AtsVj  ^tRj   eq  (16)  AtRj 

1  -0.000  056  554  796  2  -0.000  056  554  784  -0.000  056  554  785 

2  0.000  493  792  700  6  0.000  493  792  194  0.000  493  792  194 

3  0.000  259  725  929  9  0.000  259  725  401  0.000  259  725  400 

4  0.000  305  714  376  7  0.000  305  713  988  0.000  305  713  988 


46 

♦  U.S.  GOVERNMENT  PRINTING  OFFICE   I  999-0-773-022/40  1  3  I 


NIST 


Technical  Publications 

Periodical 

Journal  of  Research  of  the  National  Institute  of  Standards  and  Technology— Reports  NIST  research 
and  development  in  those  disciplines  of  the  physical  and  engineering  sciences  in  which  the  Institute  is 
active.    These  include  physics,  chemistry,  engineering,  mathematics,  and  computer  sciences.    Papers 
cover  a  broad  range  of  subjects,  with  major  emphasis  on  measurement  methodology  and  the  basic 
technology  underlying  standardization.    Also  included  from  time  to  time  are  survey  articles  on  topics 
closely  related  to  the  Institute's  technical  and  scientific  programs.   Issued  six  times  a  year. 

Nonperiodicals 

Monographs— Major  contributions  to  the  technical  literature  on  various  subjects  related  to  the  Institute's 
scientific  and  technical  activities. 

Handbooks— Recommended  codes  of  engineering  and  industrial  practice  (including  safety  codes) 
developed  in  cooperation  with  interested  industries,  professional  organizations,  and  regulatory  bodies. 

Special  Publications— Include  proceedings  of  conferences  sponsored  by  NIST,  NIST  annual  reports,  and 
other  special  publications  appropriate  to  this  grouping  such  as  wall  charts,  pocket  cards,  and 
bibliographies. 

Applied  Mathematics  Series— Mathematical  tables,  manuals,  and  studies  of  special  interest  to  physicists, 
engineers,  chemists,  biologists,  mathematicians,  computer  programmers,  and  others  engaged  in  scientific 
and  technical  work. 

National  Standard  Reference  Data  Series— Provides  quantitative  data  on  the  physical  and  chemical 
properties  of  materials,  compiled  from  the  world's  literature  and  critically  evaluated.    Developed  under  a 
worldwide  program  coordinated  by  NIST  under  the  authority  of  the  National  Standard  Data  Act  (Public 
Law  90-396).   NOTE:   The  Journal  of  Physical  and  Chemical  Reference  Data  (JPCRD)  is  published  bi- 
monthly for  NIST  by  the  American  Chemical  Society  (ACS)  and  the  American  Institute  of  Physics  (AIP). 
Subscriptions,  reprints,  and  supplements  are  available  from  ACS,  1155  Sixteenth  St.,  NW,  Washington, 
DC  20056. 

Building  Science  Series— Disseminates  technical  information  developed  at  the  Institute  on  building 
materials,  components,  systems,  and  whole  structures.   The  series  presents  research  results,  test 
methods,  and  performance  criteria  related  to  the  structural  and  environmental  functions  and  the  durability 
and  safety  characteristics  of  building  elements  and  systems. 

Technical  Notes— Studies  or  reports  which  are  complete  in  themselves  but  restrictive  in  their  treatment 
of  a  subject.   Analogous  to  monographs  but  not  so  comprehensive  in  scope  or  definitive  in  treatment  of 
the  subject  area.    Often  serve  as  a  vehicle  for  final  reports  of  work  performed  at  NIST  under  the 
sponsorship  of  other  government  agencies. 

Voluntary  Product  Standards— Developed  under  procedures  published  by  the  Department  of  Commerce 
in  Part  10,  Title  15,  of  the  Code  of  Federal  Regulations.   The  standards  establish  nationally  recognized 
requirements  for  products,  and  provide  all  concerned  interests  with  a  basis  for  common  understanding  of 
the  characteristics  of  the  products.    NIST  administers  this  program  in  support  of  the  efforts  of  private- 
sector  standardizing  organizations. 

Consumer  Information  Series— Practical  information,  based  on  NIST  research  and  experience,  covering 

areas  of  interest  to  the  consumer.    Easily  understandable  language  and  illustrations  provide  useful 

background  knowledge  for  shopping  in  today's  technological  marketplace. 

Order  the  above  NIST  publications  from:  Superintendent  of  Documents,  Government  Printing  Office, 

Washington,  DC  20402. 

Order  the  following  NIST  publications— FIPS  and  NISTIRs—from  the  National  Technical  Information 

Service,  Springfield,  VA  22161. 

Federal  Information  Processing  Standards  Publications  (FIPS  PUB)— Publications  in  this  series 
collectively  constitute  the  Federal  Information  Processing  Standards  Register.    The  Register  serves  as  the 
official  source  of  information  in  the  Federal  Government  regarding  standards  issued  by  NIST  pursuant  to 
the  Federal  Property  and  Administrative  Services  Act  of  1949  as  amended,  Public  Law  89-306  (79  Stat. 
1127),  and  as  implemented  by  Executive  Order  11717  (38  FR  12315,  dated  May  11,  1973)  and  Part  6  of 
Title  15  CFR  (Code  of  Federal  Regulations). 

NIST  Interagency  Reports  (NISTIR)— A  special  series  of  interim  or  final  reports  on  work  performed  by 
NIST  for  outside  sponsors  (both  government  and  non-government).    In  general,  initial  distribution  is 
handled  by  the  sponsor;  public  distribution  is  by  the  National  Technical  Information  Service,  Springfield, 
VA  22161,  in  paper  copy  or  microfiche  form. 


U.S.  Department  of  Commerce 

National  Institute  of  Standards  and  Technology 

325  Broadway 

Boulder,  Colorado  80303-3328 


Official  Business 

Penalty  for  Private  Use,  $300 


